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Abstract 

We  model  the  orbital  debris  environment  by  a  set  of  differential  equations  with  parameter  values  that  capture  many  of  the  complex¬ 
ities  of  existing  three-dimensional  simulation  models.  We  compute  the  probability  that  a  spacecraft  gets  destroyed  in  a  collision  during  its 
operational  lifetime,  and  then  define  the  sustainable  risk  level  as  the  maximum  of  this  probability  over  all  future  time.  Focusing  on  the 
900-  to  1000-km  altitude  region,  which  is  the  most  congested  portion  of  low  Earth  orbit,  we  find  that  -  despite  the  initial  rise  in  the  level 
of  fragments  -  the  sustainable  risk  remains  below  10-3  if  there  is  high  (>98%)  compliance  to  the  existing  25-year  postmission  deorbiting 
guideline.  We  quantify  the  damage  (via  the  number  of  future  destroyed  operational  spacecraft)  generated  by  past  and  future  space  activ¬ 
ities.  We  estimate  that  the  2007  FengYun  1C  antisatellite  weapon  test  represents  «  1%  of  the  legacy  damage  due  to  space  objects  having  a 
characteristic  size  of  >  10  cm,  and  causes  the  same  damage  as  failing  to  deorbit  2.6  spacecraft  after  their  operational  life.  Although  the 
political  and  economic  issues  are  daunting,  these  damage  estimates  can  be  used  to  help  determine  one-time  legacy  fees  and  fees  on  future 
activities  (including  deorbit  noncompliance),  which  can  deter  future  debris  generation,  compensate  operational  spacecraft  that  are 
destroyed  in  future  collisions,  and  partially  fund  research  and  development  into  space  debris  mitigation  technologies.  Our  results  need 
to  be  confirmed  with  a  high-fidelity  three-dimensional  model  before  they  can  provide  the  basis  for  any  major  decisions  made  by  the  space 
community. 

©  2009  COSPAR.  Published  by  Elsevier  Ltd.  All  rights  reserved. 
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1.  Introduction 

Orbital  debris  generated  by  50  years  of  space  activities 
poses  a  risk  for  operational  spacecraft,  which  can  collide 
in  a  catastrophic  manner  with  either  another  large  object 
(e.g.,  an  upper  stage  rocket  body)  or  with  a  smaller  frag¬ 
ment  generated  by  a  previous  collision  or  by  a  previous 
explosion  of  a  large  object  (Liou  and  Johnson,  2008).  A 
high-fidelity  three-dimensional  simulation  model  of  low 
Earth  orbit  (LEO,  which  is  the  region  between  200  and 
2000  km  altitude)  predicts  that  -  even  with  no  future 
launches  -  the  growth  rate  of  collisional  debris  would 
exceed  the  natural  decay  rate  in  «50  years  (Liou  and  Johnson, 
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2008).  Moreover,  the  analysis  in  (Liou  and  Johnson,  2008) 
did  not  account  for  the  Chinese  anti-satelhte  weapon 
(ASAT)  test  that  destroyed  the  FengYun  1C  spacecraft  in 
January  2007,  which  created  the  largest  manmade  orbital 
debris  cloud  in  history  (Liou  and  Portman,  2007).  NASA’s 
safety  guidelines  recommend  limiting  the  postmission  life¬ 
time  of  spacecraft  or  upper  stages  in  LEO  to  25  years 
(NASA  Safety  Standard  1740.14,  1995).  Because  this  mea¬ 
sure  will  not  prevent  a  positive  growth-rate  of  debris  (Liou 
and  Johnson,  2005),  it  has  been  suggested  that  the  removal 
of  large  intact  satellites  from  space  is  also  necessary  (Liou 
and  Johnson,  2006).  Although  the  impact  of  satellite 
removal  has  been  assessed  (Liou  and  Johnson,  2007),  cur¬ 
rently  there  are  no  technologies  that  are  technically  feasible 
and  economically  viable  (Liou  and  Johnson,  2006). 

Space  debris  represents  a  textbook  example  of  environ¬ 
mental  economics  (Perman  et  al.,  2003):  space  is  a  public 
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good  (i.e.,  despite  the  1976  Bogota  Declaration,  in  which 
eight  equatorial  countries  claimed  sovereignty  over  the 
portion  of  geosynchronous  Earth  orbit  lying  above  their 
territory  (Soroos,  1982),  there  are  no  well-defined  and 
enforceable  property  rights  and  no  countries  are  excluded 
from  launching  satellites)  and  debris  is  a  pollutant.  More 
specifically,  LEO  is  a  renewable  stock  resource  (much  like 
air  or  water),  in  that  debris  eventually  dissipates,  albeit  on 
an  extremely  slow  time  scale. 

The  two  fundamental  issues  in  environmental  economics 
are  to  determine  the  target  level  of  pollution  and  to  decide 
how  to  achieve  this  target.  Regarding  the  first  issue,  we  use 
the  concept  of  sustainability  (loosely  defined  as  the  highest 
utility  that  can  be  maintained  for  all  future  time),  which 
has  gained  some  popularity  in  recent  years  as  an  alternative 
to  economic  efficiency  (Perman  et  al.,  2003,  Chapter  4)  par¬ 
ticularly  for  studying  resources  that  -  like  outer  space  - 
have  no  substitutes  and  are  in  the  infancy  of  their  exploita¬ 
tion.  More  specifically,  we  introduce  the  maximum  (over 
all  future  time)  lifetime  risk  to  an  operational  spacecraft 
(i.e.,  max,>0^°(t)  in  Eq.  (14))  as  a  measure  of  sustainability. 

Regarding  the  second  fundamental  issue,  there  are  a 
number  of  different  instruments  used  in  environmental  eco¬ 
nomics  to  achieve  the  target  pollution  level:  technology 
controls,  ceilings  or  taxes  on  emissions,  subsidies  for  pollu¬ 
tion  reduction,  tradeable  emissions  permits,  and  non-com¬ 
pliance  fees.  In  practice,  technology  controls  are  often  used 
because  of  the  difficulty  of  measuring  emissions,  pollution 
levels,  or  costs,  and  indeed  the  25-year  postmission  deorbit 
guideline  and  the  end-of-life  passivation  of  rocket  bodies 
and  spacecraft  (NASA  Safety  Standard  1740.14,  1995) 
are  forms  of  technology  control. 

In  this  paper,  we  develop  an  ordinary  differential  equa¬ 
tion  (ODE)  model  that  is  a  mean  field  approximation,  incor¬ 
porating  physical  principles  for  collision  rates  (Alberty  and 
Silbey,  1997),  decay  rates  (Rossi  et  al.,  1994),  and  orbital  tra¬ 
jectories  (King-Hele,  1964),  as  well  as  rocket  body  and 
spacecraft  characteristics  (Kessler,  2000)  and  the  fragmenta¬ 
tion  models  and  empirical  distributions  developed  in  John¬ 
son  et  al.  (2001).  The  model  elucidates  the  nature  of  the 
(>  10  cm)  debris  dynamics  over  both  the  medium  (centu¬ 
ries)  and  long  (millennia)  term,  and  allows  us  to  assess  sus¬ 
tainability.  We  also  use  the  model  to  estimate  the  damage 
(in  terms  of  the  future  operational  spacecraft  that  are 
destroyed  in  collisions)  associated  with  a  present-day  launch 
or  ASAT  explosion,  as  well  as  the  damage  due  to  existing 
space  objects.  These  damages  provide  a  basis  for  quantifying 
one-time  country-specific  legacy  costs  and  setting  fees  for 
launches  and  ASAT  tests,  which  can  reduce  the  amount  of 
debris  generation,  create  a  compensation  fund  for  destroyed 
spacecraft,  and  partially  fund  the  research  and  development 
necessary  to  remove  large  objects  from  space. 

2.  The  model 

We  study  the  collisional  interactions  among  satellites 
between  the  altitudes  900  and  1000  km,  ignoring  the  effects 


of  collisions  in  lower  or  higher  orbits;  we  refer  to  this 
(900, 1000)  km  region  of  space,  which  is  the  location  of 
nearly  60%  of  all  predicted  catastrophic  colhsions  in  LEO 
over  the  next  200  years  (Liou  and  Johnson,  2008),  as  the 
shell  of  interest  (SOI).  Our  model  tracks  the  number  of 
rocket  bodies  and  three  types  of  spacecraft  over  time, 
and  the  effective  number  of  four  types  of  fragments.  The 
effective  number  weights  each  fragment  by  the  fraction  of 
time  its  orbit  is  in  the  SOI  and  the  probability  that  it  is  haz¬ 
ardous  or  benign  (as  explained  below);  spacecraft  and 
rocket  bodies,  which  together  are  referred  to  as  intacts, 
are  assumed  to  have  circular  orbits,  and  hence  spend  all 
of  their  orbital  time  in  the  SOI.  Because  our  model  treats 
the  SOI  in  isolation,  we  perform  a  sensitivity  analysis  in 
Section  3.2  that  allows  fragments  created  in  other  shells 
to  enter  the  SOI. 

Some  intacts  in  LEO  have  a  deorbit  capability  to  facili¬ 
tate  a  25-year  postmission  lifetime  in  LEO  in  accordance 
with  NASA  guidelines  (NASA  Safety  Standard  1740.14, 
1995):  at  the  end  of  its  operational  lifetime,  such  an  intact 
immediately  maneuvers  to  either  an  orbit  with  a  lower  peri¬ 
gee  (typically  below  600  km)  from  which  it  decays  naturally 
from  LEO  within  25  years,  or  to  a  disposal  orbit  above 
2500  km  altitude  (intacts  operating  in  the  SOI  will  typically 
choose  the  former  option).  By  deorbit,  we  shall  mean  that 
the  intact  leaves  the  SOL  Consequently,  we  differentiate 
the  total  number  of  spacecraft  ( S(t ))  into  three  types:  those 
that  have  a  deorbit  capability  ( Sd{t )),  those  that  do  not  and 
are  still  operational  (S® (?)).,  and  those  that  do  not  have  a 
deorbit  capability  and  are  no  longer  operational  (S„(f)). 
(We  use  the  same  symbol  for  two  purposes:  e.g.,  Sd(t)  is 
the  effective  number  of  satellites  of  type  Sd  at  time  t.)  Future 
spacecraft  are  launched  into  the  SOI  at  annual  rate  and 
we  assume  that  some  recent  spacecraft  and  a  fraction  6d  of 
all  future  spacecraft  have  a  postmission  deorbit  capability. 
These  spacecraft  deorbit  at  rate  nr„  where  is  the  average 
operational  lifetime.  The  majority  of  existing  spacecraft  in 
LEO  and  a  fraction  1  —  9d  of  future  spacecraft  do  not  have 
a  postmission  deorbit  capability,  and  change  category,  from 
S°  to  Sn,  at  rate  firi,  and  then  deorbit  naturally  at  rate  jin. 
Similarly,  some  upper  stage  rocket  bodies  immediately 
deorbit  the  SOI  and  do  not  appear  in  our  model,  while  other 
rocket  bodies  {R(t))  linger  at  the  separation  altitude  and 
deorbit  naturally  at  rate  nR;  these  latter  rocket  bodies  are 
inserted  into  the  SOI  at  the  annual  rate  ).R.  Although  deor¬ 
biting  of  spacecraft  and  rocket  bodies  may  actually  take 
several  months,  sensitivity  analysis  (details  not  shown) 
reveals  that  our  assumption  of  instantaneous  deorbiting 
has  a  negligible  impact  on  all  of  our  results. 

We  assume  passivation  techniques  (e.g.,  venting  rocket 
bodies  of  residual  fuel  and  discharging  spacecraft  batteries) 
are  in  place  to  prevent  future  explosions,  and  so  debris  is 
generated  in  our  model  only  via  colhsions  between  satel¬ 
lites.  Catastrophic  collisions  result  in  the  fragmentation 
of  both  objects,  whereas  noncatastrophic  collisions  result 
in  the  fragmentation  of  the  less  massive  object.  Essentially 
all  intact-intact  collisions  are  catastrophic.  We  consider 
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fragments  having  a  characteristic  length  >  10  cm,  gener¬ 
ally  the  lower  limit  of  what  sensors  in  the  U.S.  Space  Sur¬ 
veillance  Network  can  resolve  (Johnson  et  al.,  2001).  We 
refer  to  fragments  as  hazardous  or  benign  depending  on 
whether  or  not  they  can  produce  catastrophic  collisions 
with  intacts;  although  fragments  can  be  hazardous  or 
benign  to  each  other  as  well,  this  distinction  is  far  less 
important  because  many  fewer  fragments  are  generated 
in  fragment-fragment  collisions.  Hence,  the  effective  frag¬ 
ment  number  is  indexed  by  two  symbols  and  will  be  written 
Fz,  where  k  e  {/?,  b]  refers  to  being  hazardous  or  benign  to 
intacts  and  %  e  {R,S}  refers  to  the  source  of  a  fragment  - 
rocket  body  or  spacecraft.  A  particular  fragment  is  not 
simply  hazardous  or  benign  to  an  intact:  the  uncertainty 
in  collision  velocity  causes  the  properties  of  the  fragment 
to  determine  the  probability  with  which  it  is  hazardous 
or  benign  in  a  particular  collision.  A  particular  fragment 
from  source  i  e  {R, S'},  therefore,  increases  the  effective 
number  of  Fb  and  Fbz  by  quantities  that  sum  to  one.  Frag¬ 
ments  deorbit  at  constant  rate  /v- 

Satellites  of  types  a  and  y  participate  in  collisions  at  rate 
/?ayoc(f)y(f)  at  time  t,  where  a,y  e  {.S’,  R,Fz}-  hence,  is 
the  number  of  collisions  between  satellites  of  types  oc  and 
y  per  unit  time  per  satellite  of  type  a  per  satellite  of  type 
y.  Because  a  satellite’s  cross  section  (which  dictates  its  like¬ 
lihood  of  collision)  and  mass  (which  determines  the  num¬ 
ber  of  new  fragments  if  it  is  involved  in  a  collision)  have 
a  joint  probability  distribution,  we  cannot  estimate  in  iso¬ 
lation  the  number  of  fragments  resulting  from  a  collision 
between  two  types  of  satellites.  Rather,  we  assume  that 
fragments  of  type  F *  are  generated  from  collisions  between 
satellites  of  types  a  and  y  at  rate  8™x(t)y(t)  at  time  t,  for 
a,  y  6  {5,  R,  Fz}  (a  factor  \  is  necessary  if  a  =  y).  For  exam¬ 
ple,  a  collision  between  a  rocket  body  and  a  spacecraft  pro¬ 
duces  all  four  types  of  fragments:  hazardous  and  benign 
fragments  from  both  a  rocket  body  and  a  spacecraft. 

The  effective  number  of  spacecraft,  hazardous  fragments, 
benign  fragments,  and  total  fragments  at  time  t  are  given  by, 
respectively,  S(f)  =  S°(t)  +  Sn(t)+  Sd(t),Fh(t)  =  FR(t)  +Fhs 
(t),Fb(t)=FbR(t)+Fb(t),  and  Fs(t)=  Fh{i)+Fb{t).  Let 
U  =  {S,R,Fb,Fb,FhR,FhR\  be  the  set  of  satellite  types, 
Uh  =  {5,  R,  Fb.  FR\  be  the  set  of  satellite  types  hazardous 
to  intacts,  UF  =  {Fb.  Fb,  FR,  FbR)  be  the  set  of  fragment 
types,  and  U1  =  {.S’,  R}  be  the  set  of  intact  types.  The  model 
is  given  by 


*(0  =  ^-£M0«M-MW,  (1) 

a  eUh 

s°n(t)  =  (1  -  ed)xB  -  £  PsJWWt)  -  HoS°n(t),  (2) 

a  eU* 

S.(t)  =  HoS°n(t)  -  £  PsMiHt)  -  finSn(t),  (3) 

a  £Uh 

Sd(t)  =  OAo  -  £  -  nMt ),  (4) 

<xeuh 

*?(<)  =  \  £  £ (5) 

Z  a<=E/  yeU 


for  t  6  U1,  k  e  {h,  b}. 

A  key  shortcoming  of  the  model  relates  to  our  handling 
of  fragment-fragment  collisions.  We  assume  that  a  frag¬ 
ment  generated  in  a  fragment-fragment  collision  inherits 
the  same  properties  as  its  ancestral  intact  (i.e.,  rocket  body 
or  spacecraft).  Because  a  fragment  generated  by  a  frag¬ 
ment-fragment  collision  would  be  smaller  than  its  parent 
fragments,  we  view  Eq.  (5)  as  upper  bounding  the  effect 
of  fragment-fragment  collisions,  and  then  as  a  lower 
bound  consider  a  variant  of  Eq.  (5)  that  simply  disregards 
fragment-fragment  collisions  (i.e.,  sets  S™  =  0  for 
a,  y  £  UF).  We  will  usually  be  ignoring  fragment-fragment 
collisions.  Hence,  we  refer  to  the  version  of  the  model  that 
incorporates  fragment-fragment  collisions  (i.e.,  with 
(5“  >0  for  oi,  y  €  UF)  as  the  fragment-fragment  version 
of  the  model;  otherwise,  the  reader  should  assume  that 
we  are  ignoring  fragment-fragment  collisions  (i.e.,  setting 
d “  =  0  for  a,  y  6  UF).  The  intact-intact  parameter  values 
are  the  same  in  both  versions  of  the  model.  When  we  ignore 
fragment-fragment  collisions,  the  intact-benign  fragment 
parameter  values  account  only  for  the  loss  of  a  benign  frag¬ 
ment,  and  the  intact-hazardous  fragment  parameter  values 
account  only  for  the  breakup  of  the  intact  and  the  loss  of  a 
hazardous  fragment. 

2.1.  Equilibrium  analysis 

The  main  goal  of  this  subsection  is  to  understand  the 
long-term,  or  equilibrium,  behavior  of  system  (l)-(5). 
Toward  this  end,  we  consider  a  simpler  model  that  tracks 
only  two  types  of  satellites,  intacts  (7(f))  and  hazardous  frag¬ 
ments  {F(t)),  and  ignores  fragment-fragment  collisions: 


'm = -pni2{i)  -  PiFmm  -  + xh 

(6) 

Ht)  =  <V2M  +  <W(f)F(t)  -  iiFF{t). 

(7) 

We  define  the  parameters  in  Eqs.  (6)  and  (7) 
those  in  Eqs.  (1)— (5)  via 

in  terms  of 

=  max  flxy, 

(8) 

PIF  =  max 

(9) 

Hj  =  min{{if  XR  >  0,  nR,  else  oo}, 

(10) 

{if  9d  <  1,  jU„;  else  oo },n0}, 

Xi  =  XR  +  X0, 

(11) 

bn  8if  f  d  C ) 

s=&=manmaxwmaxs)’ 

(12) 

HF  =  min^, 

for  a, ye  U1,  C,qeUF,Te  {R,  S},  k  e  {h,  b}. 

(13) 

In  the  remainder  of  this  subsection,  we  state  mathemat¬ 
ical  results  that  are  rigorously  derived  in  Section  B  of  Brad¬ 
ley  and  Wein  (2008),  which  is  an  earlier  version  of  this 
paper. 

First,  if  7(0)  =  T’(O),  where  P'(t)  =  R(t)  +  S°(t)  + 
Sn(t )  +  Sd(t)  is  the  total  number  of  intacts  in  Eqs.  ( 1 )— (5), 
and  F(0)  =  Fs(0),  then  F(t)  >  Fs(t)  for  t  >  0;  i.e.,  the 


A.M.  Bradley,  L.M.  Weinl Advances  in  Space  Research  43  (2009)  1372-1390 


1375 


number  of  hazardous  fragments  in  Eqs.  (6)  and  (7)  provides 
an  upper  bound  on  the  number  of  fragments  in  Eqs.  (1)- 

(5) . 

Second,  the  number  of  hazardous  fragments  F(t)  in  Eqs. 

(6)  and  (7)  remains  bounded,  and  hence  the  total  number  of 
fragments  Fs(t)  in  (1)— (5)  also  remains  bounded;  this  is  the 
main  result  of  this  subsection.  We  also  derive  explicit 
expressions  for  the  equilibrium  number  of  intacts  and  the 
equilibrium  number  of  fragments  for  Eqs.  (6)  and  (7). 

Third,  if  all  space  participants  fully  comply  with  deorbit 
procedures  for  rocket  bodies  and  spacecraft  (i.e.,  kR  =  0 
and  0d  =  1),  Eqs.  (l)-(7)  are  identical  as  t  — >  oo  with  iden¬ 
tical  parameters  (e.g.,  d„.  =  dfFh),  with  the  exception  that 
the  factor  \  is  absorbed  into  the  parameter  bu. 

Fourth,  if  we  include  the  fragment-fragment  term 
SffF 2(t)  in  Eq.  (7),  then  the  upper  bounding  of  Fs(t)  in 
Eqs.  (1)— (5)  by  F{t)  in  Eqs.  (6)  and  (7)  still  holds.  There 
are  two  possibilities  as  t  — >  oo:  the  number  of  hazardous 
fragments  grows  without  bound,  or  the  model  approaches 
an  equilibrium  that  is  nearly  identical  to  that  of  the  model 
without  fragment-fragment  collisions  (i.e.,  in  which 
bFF  =  0)  (this  case  occurs  with  our  set  of  parameter  values 
under  full  compliance). 

2.2.  Parameter  estimation 

Values  for  the  parameters  in  Eqs.  (1)— (5)  are  estimated 
in  (and  displayed  in)  the  Appendix  using  Monte  Carlo  inte¬ 
gration.  Using  the  empirical  probability  distributions  in 
Johnson  et  al.  (2001),  we  assign  each  fragment  a  random 
characteristic  length,  area-to-mass  ratio,  and  magnitude 
and  direction  of  change  in  velocity  after  a  collision;  the  first 
two  random  variables  dictate  a  fragment’s  area  and  mass, 
while  the  third  combines  with  the  first  two  to  give  the 
hazard  probability.  In  addition,  each  fragment’s  source’s 
altitude  is  uniformly  distributed  within  the  SOI,  and  the 
direction  of  its  post-collision  change  in  velocity  is  uni¬ 
formly  distributed  on  the  surface  of  the  sphere.  The  colli¬ 
sion  rate  parameters  {Pxy),  the  debris  generation 
parameters  (c>™),  and  the  decay  rates  are  computed  by  tak¬ 
ing  expectations  over  the  three  probability  distributions 
from  Johnson  et  al.  (2001)  and  the  two  uniform 
distributions. 

In  the  simulation,  each  fragment  is  weighted  by  the  frac¬ 
tion  of  time  it  resides  in  the  SOI,  which  is  derived  using  the 
equations  of  orbital  mechanics  (e.g.,  King-Hele,  1964, 
Chapter  3).  In  addition,  a  fragment  from  source 
t  e  {R.  5}  is  split  between  hazardous  (Fhr)  and  benign 
(Fb)  by  calculating  the  probability  that  a  collision  between 
it  and  an  intact  would  be  catastrophic,  which  occurs  if  the 
kinetic  energy  of  the  projectile  in  Joules  is  40  times  greater 
than  the  mass  of  the  target  in  grams  (Johnson  et  al.,  2001). 

To  estimate  the  collision  rate  parameters,  we  use  the 
ideal  gas  model  (Alberty  and  Silbey,  1997)  to  derive  an 
analytical  expression  for  in  terms  of  the  mean  relative 
velocity  of  satellites  in  the  SOI  (Johnson  et  al.,  2001),  the 
volume  of  the  SOI,  and  the  mean  collision  cross  sections 


of  the  various  satellite  types  (using  data  from  Kessler 
(2000)  for  intacts  and  the  above  probability  distributions 
for  fragments).  We  apply  the  ideal  gas  model  to  each  por¬ 
tion  of  the  SOI  at  a  particular  latitude,  and  use  the  full  cat¬ 
alog  of  satellite  element  sets  from  Space  Track  (2008)  to 
derive  an  empirical  distribution  of  effective  spatial  density 
of  satellites  vs.  latitude  to  capture  the  fact  that  the  spatial 
density  of  satellites  -  and  hence  the  collision  rate  -  is  higher 
near  the  poles. 

To  estimate  the  debris  generation  parameters,  we  use 
Eq.  (4)  in  Johnson  et  al.  (2001),  which  specifies  the  number 
of  fragments  from  a  collision  having  characteristic  length 
of  at  least  a  certain  size,  in  terms  of  the  mass  of  the  two 
objects  involved  in  the  collision,  the  collision  velocity, 
and  whether  the  collision  was  catastrophic. 

We  estimate  the  decay  rates  as  reciprocals  of  residence 
times.  We  set  p~l  =  3  years  to  reflect  the  historical  average 
of  a  3-year  mission  lifetime  and  the  25-year  decay  rule  rec¬ 
ommended  by  NASA  safety  guidelines  (NASA  Safety  Stan¬ 
dard  1740.14,  1995,  2005).  To  estimate  pR,p„,  and  pF* ,  we 
use  data  in  Table  4  of  Rossi  et  al.  (1994)  to  equate  the  res¬ 
idence  time  of  a  satellite  to  1 10  years  divided  by  the  average 
area-to-mass  ratio  for  the  satellite  type,  which  is  quantified 
using  the  probability  distributions  in  Johnson  et  al.  (2001) 
for  fragments  and  the  values  in  Kaula  (1983)  for  intacts. 

To  estimate  the  launch  rate  of  new  spacecraft  into  the 
SOI,  we  analyze  the  full  catalog  of  satellite  element  sets 
at  the  end  of  2007  from  Space  Track  (2008),  identify  the 
satellites  that  were  spacecraft,  and  determine  the  effective 
number  of  spacecraft  in  the  SOI.  Inspecting  the  annual 
launch  rates  of  these  spacecraft  into  the  SOI  over  the  past 
decade,  we  set  ka  =  3/year.  A  similar  analysis  for  rocket 
insertions  leads  to  kR  =  1/year;  because  only  a  few  rocket 
launchers  remain  that  do  not  enable  rocket  body  deorbit¬ 
ing,  we  assume  kR  =  1/year  for  10  years  and  kR  =  0  thereaf¬ 
ter.  While  the  spacecraft  deorbiting  guidelines  have  been 
implemented,  they  have  not  been  fully  adopted  by  major 
space  agencies  (Liou  and  Johnson,  2007),  and  we  set 
0d  =  |  and  perform  a  sensitivity  analysis  on  6d  later. 

The  full  catalog  of  satellite  element  sets  from  Space 
Track  (2008)  are  also  used  to  determine  the  initial  condi¬ 
tions  for  rocket  bodies,  spacecraft,  and  fragments: 
i?(0)  =  183.3,5(0)  =207.2,  and  Fs(0)  =  955.5.  We  set 
S"(0)  =  =  3  and  Sd(0)  =  ^  =  6,  which  are  the  equi¬ 

librium  values  if  we  ignore  collisions  in  Eqs.  (2)  and  (4),  and 
S„(0)  =  S(0)  -  S°(0)  -  Sd(0).  To  allocate  Fs(0)  among  the 
four  fragment  types,  we  assume  that  -  aside  from  the 
FengYun  1C  fragments,  which  are  known  to  originate  from 
a  spacecraft  -  90%  of  the  initial  fragments  are  from  rocket 
bodies  (our  results  are  insensitive  to  this  assumption)  and 
use  the  simulated  proportion  of  fragments  of  type  Fz  gener¬ 
ated  by  the  breakup  of  an  intact  (Section  F). 

2.3.  Comparison  to  Liou  and  Johnson  (2008) 

NASA  researchers  have  recently  developed  a  three- 
dimensional  simulation  of  near-Earth  orbital  space  that 
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tracks  each  satellite  over  time  (Johnson  et  al.,  2001;  Liou 
and  Johnson,  2005,  2008).  Satellite  breakups  produce  frag¬ 
ments  whose  properties  are  drawn  from  the  distributions 
described  in  Johnson  et  al.  (2001),  though  with  some  mod¬ 
ifications  as  described  in  Liou  and  Johnson  (2008).  We 
compare  the  results  of  the  fragment-fragment  version  of 
our  model  against  those  presented  in  Liou  and  Johnson 
(2008)  for  the  SOI  (Nicholas  Johnson,  private  communica¬ 
tion),  which  provide  mean  values  for  150-  to  200-year  sim¬ 
ulation  runs.  We  set  =  1*  =  S°(0)  =  Srf(0)  =  0  to 
coincide  with  the  assumptions  in  Liou  and  Johnson 
(2008):  no  future  launches,  full  compliance  with  rocket 
body  deorbiting,  and  all  initial  spacecraft  are  nonopera- 
tional  and  have  no  deorbit  capability.  The  initial  condi¬ 
tions  evident  in  Fig.  8  of  Liou  and  Johnson  (2008)  are 
«  490  “intacts  +  mission-related  debris”  and  520  frag¬ 
ments.  It  is  not  indicated  what  proportion  of  intacts  are 
rocket  bodies.  Because  intact  rocket  bodies  are  larger  and 
more  spatially  clustered  (Section  C.3)  than  intact  space¬ 
craft,  we  must  take  care  in  determining  our  initial  condi¬ 
tions  for  the  comparison  with  Liou  and  Johnson  (2008). 
In  our  initial  conditions  for  the  comparison,  we  set  the  ini¬ 
tial  number  of  rocket  bodies  to  183.3  (Space  Track,  2008), 
the  initial  number  of  spacecraft  to  490  -  183.3  =  306.7, 
and  the  initial  number  of  fragments  to  520,  which  are  allo- 


Year 


Fig.  1.  Comparison  of  the  fragment-fragment  version  of  Eqs.  (1)— (5)  ( — ) 

to  the  results  in  Fig.  8  of  Liou  and  Johnson  (2008)  ( - ).  The  model  of 

Liou  and  Johnson  (2008)  does  not  distinguish  between  hazardous  and 
benign  fragments. 


cated  among  the  four  categories  using  the  methods 
described  earlier. 

Fragment  growth  in  both  models  is  nearly  linear  for  the 
first  200  years  (Fig.  1),  and  our  model  underestimates  the 
number  of  fragments  after  200  years  relative  to  Liou  and 
Johnson  (2008)  by  8.6%  (2356  vs.  2579).  Our  model  overes¬ 
timates  the  number  of  destroyed  intacts  (1 1.02  vs.  8.21)  and 
the  number  of  intact-intact  collisions  (4.38  vs.  2.65),  under¬ 
estimates  the  fraction  of  intact-fragment  collisions  that  are 
catastrophic  (0.33  vs.  0.44),  and  overestimates  the  fraction 
of  fragment-fragment  collisions  that  are  catastrophic  (0.81 
vs.  0.67)  (Table  1).  Possible  reasons  for  these  discrepancies 
are  given  in  Section  4. 

3.  Results 

3.1.  Base-case  results 

Numerical  solutions  to  both  versions  of  Eqs.  (1)— (5) 
(Fig.  2)  reveal  that  in  both  cases,  the  total  number  of  frag¬ 
ments  increases  in  a  convex  manner  for  «  1500  years.  In 
the  absence  of  fragment-fragment  collisions,  the  total 
number  of  fragments  approaches  its  asymptotic  value 
[Fs  =  4.7  x  1 05 ,  F*  =  3.1  x  105),  while  in  the  presence  of 
fragment-fragment  collisions  the  number  of  fragments 
grows  without  bound,  blowing  up  at  time  («  1473 
(although  in  another  set  of  simulations  we  found  that  the 
fragment  level  remains  bounded  if  all  spacecraft  deorbit 
after  their  mission,  i.e.,  if  6d  =  1).  The  two  versions  of 
the  model  are  nearly  indistinguishable  for  the  first 
500  years  (e.g.,  at  200  years,  7^(200)  =  1006  and  1015  in 
the  absence  and  presence  of  fragment-fragment  collisions). 
Hereafter,  we  assume  there  are  no  fragment-fragment 
collisions. 

The  catastrophic  intact-fragment  collisions  are  the  least 
frequent  type  of  collision  over  the  first  several  hundred 
years,  but  eventually  become  the  most  common  type  of  col¬ 
lision  (Fig.  3a  and  b).  Our  primary  metric  for  measuring 
the  hazard  is  the  probability  that  a  spacecraft  launched 
at  time  t  will  be  destroyed  (via  an  intact-intact  or  cata¬ 
strophic  intact-fragment  collision)  while  it  is  still  opera¬ 
tional,  which  is 

^(0  =  1- (i -Emo|  •  (i4) 


Table  1 

Comparison  of  the  fragment-fragment  version  of  Eqs.  (1)— (5)  to  Liou  and  Johnson  (2008).  Lor  Liou  and  Johnson  (2008),  these  are  the  average  number  of 
collisions  (averaged  over  150  Monte  Carlo  simulations)  in  the  SOI  for  the  results  appearing  in  Fig.  8  in  Liou  and  Johnson  (2008)  (Nicholas  Johnson, 
private  communication). 


Object  type 

Number  of  collisions  for  the  next  200  years 

Catastrophic/  non-catastrophic 

(1M5) 

Liou  and  Johnson  (2008) 

(1H5) 

Liou  and  Johnson  (2008) 

Intact-intact 

4.38 

2.65 

4.38/0.00 

2.64/0.01 

Intact-fragment 

6.83 

6.43 

2.26/4.57 

2.83/3.60 

F  ragment-fragment 

0.58 

0.58 

0.47/0.11 

0.39/0.19 

Total 

11.79 

9.66 

7.11/4.68 

5.06/3.27 
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Fig.  2.  The  solution  to  both  versions  of  Eqs.  (1)— (5)  for  (a)  200  years  and 
(b)  104  years.  The  dashed  curves  ignore  fragment-fragment  collisions,  and 
the  solid  curves  include  fragment-fragment  collisions. 


Fig.  3.  (a)  Cumulative  number  of  intact-intact,  catastrophic  intact- 
fragment,  and  noncatastrophic  intact-fragment  collisions  over  the  first 

200  years  and  (b)  cumulative  number  ( — )  and  instantaneous  rate  ( - )  of 

these  collisions  over  104  years. 


Eq.  (14)  is  constructed  by  considering  the  probability  a  par¬ 
ticular  operational  spacecraft  avoids  a  collision  per  year, 
raising  the  probability  to  a  power  equal  to  the  mean  oper¬ 
ational  lifetime,  and  subtracting  the  result  from  1  to  obtain 
the  risk.  This  lifetime  risk  is  initially  1.84  x  10“4  and  in¬ 
creases  to  2.82  x  10-4  in  200  years  (Fig.  4a).  However,  this 
risk  rapidly  increases  in  sa  1000  years  for  sa  2000  years, 
eventually  approaching  the  equilibrium  value  of 
2.19  x  10‘2  (Fig.  4b). 

3.2.  Sensitivity  analyses 

Because  rocket  bodies  that  do  not  immediately  deorbit 
are  only  launched  during  the  first  10  years  in  our  model, 
our  results  are  highly  insensitive  to  the  value  of  XR.  We  vary 
the  spacecraft  deorbit  compliance  rate  9d  while  keeping  the 
other  parameters  fixed  at  their  base-case  value  (Fig.  5).  We 
define  the  sustainable  lifetime  risk  to  be  the  maximum  of 
the  lifetime  risk  over  all  future  times,  which  is  given  by 
max(>0T°(/),  where  r°(t)  is  defined  in  Eq.  (14).  While  the 
lifetime  risk  to  an  operational  spacecraft  over  the  next 
200  years  remains  below  4.2  x  10  4  regardless  of  the  value 
of  9d,  the  sustainable  lifetime  risk  decreases  to  10-2  and 


10-3,  respectively,  if  the  deorbit  compliance  increases  to 
84.9%  and  98.2%,  respectively.  The  sustainable  lifetime  risk 
drops  sharply  for  very  high  compliance  rates  and  is 
4.9  x  1 0  4  under  full  compliance  (i.e.,  9d  =  1). 

In  the  full  compliance  case,  the  fragment  level  peaks  at 
10,800,  which  is  11.3-fold  higher  than  the  current  level, 
after  several  millennia,  and  then  approaches  a  very  small 
equilibrium  value  (.P  =  2.61).  Under  full  compliance,  if 
the  spacecraft  launch  rate  is  increased  10-fold  (e.g.,  due 
to  the  emergence  of  reusable  launch  vehicles),  then  the  life¬ 
time  risk  in  200  years  is  2.53  x  10-4  and  the  sustainable 
lifetime  risk  is  9.05  x  10-4. 

Our  model  treats  the  SOI  as  isolated  from  all  other  alti¬ 
tudes.  However,  fragments  often  have  eccentric  orbits  or 
decay  from  higher  altitudes,  and  so  fragments  created  out¬ 
side  the  SOI  can  interact  with  satellites  inside  it.  To  under¬ 
stand  the  effect  of  this  issue,  we  introduce  a  parameter  s 
into  (5): 

Kit)  =  E{\^Z  _ ^K{t). 

\Z  *€U  y€U  ) 


(15) 
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Fig.  4.  The  lifetime  risk  to  an  operational  spacecraft  (see  Eq.  (14)),  which  is 
the  probability  a  spacecraft  launched  at  a  particular  time  will  be  destroyed 
(via  an  intact-intact  or  a  catastrophic  intact-fragment  collision)  while  it  is 
still  operational,  plotted  vs.  time  for  (a)  200  years  and  (b)  104  years. 


Fig.  5.  The  sustainable  lifetime  risk  to  an  operational  spacecraft,  i.e. 

max(;jo r°{t)  in  Eq.  (14),  ( — )  and  the  lifetime  risk  at  t  =  200  years  ( - )  vs. 

the  proportion  6d  of  new  launches  that  deorbit  the  spacecraft  at  the  end  of  its 
operational  life. 

If  s  >  1,  then  the  number  of  fragments  created  in  collisions 
is  greater  than  in  the  nominal  case  e  =  1 .  The  new  frag¬ 
ments  may  be  thought  of  as  coming  from  a  different  shell 


having  dynamics  similar  to  those  of  the  SOI.  We  find  that 
the  sustainable  lifetime  risk  and  the  lifetime  risk  after 
200  years  are  approximately  linear  functions  of  e.  For 
e  =  2,  the  sustainable  lifetime  risk  is  2.1  times  greater, 
and  the  risk  after  200  years  is  1.2  times  greater,  than  in 
the  nominal  case.  Based  on  calculations  similar  to  those 
in  Section  B.l,  s  =  2  is  an  overestimate  of  the  value  needed 
to  account  for  fragments  from  other  shells. 

3.3.  On  the  optimality  of  full  deorbit  compliance 

Motivated  by  the  results  of  the  sensitivity  analyses,  we 
compute  the  optimal  compliance  rate  from  a  societal,  sus¬ 
tainable  perspective.  Let  Cd  be  the  cost  to  deorbit,  Cs  be 
the  cost  of  a  destroyed  operational  spacecraft,  9a  be  the 
attempted  compliance  rate,  s  be  the  probability  that  a  deor¬ 
bit  attempt  is  successful  (so  that  9d  =  s9a),  and  be 

the  sustainable  lifetime  risk  of  an  operational  spacecraft 
(i.e.,  fmJ9d)  =  max^oF'M)  when  the  successful  compli¬ 
ance  rate  is  9d.  Then  to  minimize  the  sustainable  costs 
related  to  spacecraft  deorbiting  and  destruction,  we  solve 
min  ga9aCd  +  r^mx(s9a)Cs. 

We  estimate  the  velocity  requirements  for  deorbiting 
from  the  SOI  to  be  Av  =  128  m/s  (Section  B.l).  Extrapo¬ 
lating  from  the  cost  estimates  in  Table  3  of  Wiedemann 
et  al.  (2004)  (which  are  in  FY  2002  dollars),  the  additional 
launch  cost  «  $0.5M  and  the  total  development,  produc¬ 
tion,  and  launch  cost  is  about  $20M.  We  ignore  the  devel¬ 
opment  and  production  cost,  and  assume  the  postmission 
deorbit  cost  is  Cd  =  $0.5M.  The  replacement  cost  of 
building  and  launching  a  new  satellite  is  «  $0.2-0. 5B  for 
a  commercial,  scientific,  or  weather  satellite,  or 
« $0.5-1 .5B  for  a  defense  satellite  (de  Week  et  al., 
2003),  and  we  set  Cs  =  $0.5B.  The  failure  rate  of  propul¬ 
sion  systems  has  been  estimated  to  be  3.9%  (de  Week 
et  al.,  2003,  Fig.  8)  which  yields  5  =  0.961.  Substituting 
these  values  into  the  optimization  problem  and  using 
the  F/ax(0j)  function  in  Fig.  5,  we  find  that  9a=  1  is 
optimal. 

To  assess  the  robustness  of  this  result,  we  determine  how 
high  the  deorbit  cost  has  to  be  to  make  9a  =  1  no  longer 
optimal,  which  occurs  when  either  the  derivative 
of  the  total  sustainable  cost  with  respect  to  9a  equals  zero 
at  9 a  =  1,  i.e.,  Cd  =  Cs,  or  when  the  cost 

dua  |fla=1 

at  9 a  =  0  equals  the  cost  at  9a=  1,  i.e.,  Cd  = 
(Cax(0)  -  Cax(‘s))C,5-  The  breakeven  value  is 
Cd  =  $30. 0M,  which  is  60  times  larger  than  the  actual  value 
and  several  million  dollars  more  than  the  total  cost  esti¬ 
mate  in  Wiedemann  et  al.  (2004).  Although  we  omit  the 
details,  a  similar  derivation  shows  that  100%  compliance 
of  rocket  body  deorbiting  is  also  optimal. 

3.4.  Assessing  the  damage  due  to  space  activities 

A  prerequisite  for  quantifying  fees  is  to  calculate  the 
damage  caused  by  various  space  activities.  Toward  this 
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end,  we  define  the  damage  caused  by  a  space  activity  to  be 
the  total  number  of  destroyed  operational  spacecraft  gener¬ 
ated  by  this  activity.  If  we  let  S°(t)  =  S°(t)  +  Sd{t)  be  the 
number  of  operational  spacecraft  at  time  t,  then  the  num¬ 
ber  of  destroyed  operational  spacecraft  up  to  time  T  is 

nd(T)=  [T  S°(t)J2Psyy(t)  df-  (16) 

Jo  yeUh 

Let  nd  (T)  be  the  damage  up  to  time  T  in  the  current  envi¬ 
ronment,  and  nd  (  T)  be  the  damage  up  to  time  T  after  a  per¬ 
turbation  has  been  applied  at  time  0.  We  define  an 
activity’s  damage  up  to  time  T  to  be  nd(T)  —  nd{T).  The 
perturbations  are  i?(0)  — >  i?(0)  +  1  for  an  extra,  rocket 
insertion,  R{ 0)  -*•  i?(0)  -  1  and  F\{ 0)  -f  F\{ 0)  -h^16”0^ 
for  a  rocket  breakup,  and  S„ (0)  — >  S„( 0)  —  1  and  FKS( 0)  — > 
FKS{0)  +  for  a  satellite  breakup,  where  the  latter 

two  perturbations  follow  from  Eq.  (D.ll)  and  N(-)  and 
s™  are  defined  in  Eqs.  (D.3),  (D.4),  and  (D.7).  For  a  new 
extra  launch,  we  do  not  want  to  consider  the  self-inflicted 
damage  if  the  newly-launched  spacecraft  is  destroyed  while 
it  is  operational,  although  we  do  want  to  capture  the  sub¬ 
sequent  damage  caused  by  the  fragments  of  such  a  colli¬ 
sion.  Hence,  for  a  new  extra  launch  of  a  spacecraft  with 
deorbit  capability,  we  add  the  equation  Sd(t)  = 

-  M*)  to  Eqs-  (1M5),  let  srd  be  in 

the  set  Uh,  and  consider  the  perturbation  5^(0)  — >  1  (from 
5^(0)  =  0).  For  a  new  extra  launch  of  a  spacecraft  without 
deorbit  capability,  we  simply  set  5„(0)  — ►  5„(0)  +  1  because 
an  operational  spacecraft  whose  loss  does  not  count  as  an 
operational  loss  can  be  considered  inoperational  from  the 
very  start. 

The  legacy  cost  is  the  difference  between  the  damage  in 
the  current  environment  and  the  damage  in  an  environment 
in  which  there  are  no  nonoperational  hazardous  objects; 
i.e.,  we  use  the  perturbation  F*(0)  — >  0,5„(0)  — >  0,  and 
i?(0)  — >  0.  In  addition,  we  assess  the  portion  of  the  legacy 
cost  due  to  the  FengYun  1C  ASAT  test  via  the  perturbation 
F*(0)  -»•  **(Q)  -  149.1  and  Fbs{ 0)  -►  F*(0)  -  251.8  (Sec¬ 
tion  F). 

The  damage  generated  by  a  perturbation  lasts  for 
«  104  years,  and  hence  nd(T)  -  nd{T)  approaches  a  con¬ 
stant  as  T  — »  oo  (Fig.  6).  In  the  short  term,  a  breakup  pro¬ 
duces  more  damage  than  an  intact  of  the  same  type  (Fig.  6); 
in  the  long  term,  the  opposite  is  true  because  an  intact 
remains  in  orbit  longer  than  a  fragment  and  so  is  hazardous 
for  longer.  In  keeping  with  the  sustainability  philosophy 
(long  time  horizons  are  also  considered  for  the  nuclear 
waste  depository  at  Yucca  Mountain,  which  needs  to  com¬ 
ply  for  104  years  (Long  and  Ewing,  2004)),  we  consider  an 
infinite  horizon  and  refer  to  the  asymptotic  value  of 
nd  (T)  —  nd(T)  as  the  damage.  The  damage  from  launching 
a  spacecraft  with  deorbit  capability  is  2.15  x  10-5  destroyed 
spacecraft,  which  is  3000-fold  less  damaging  than  launching 
a  spacecraft  without  deorbit  capability.  The  FengYun  1C 
ASAT  test  can  be  expected  to  cause  the  same  damage  as  8 
rocket  body  insertions  or  2.6  launches  of  spacecraft  that 


T  (years) 


Fig.  6.  The  damage  (i.e.,  destroyed  operational  spacecraft)  up  until  time  T 
caused  by  various  activities  that  occur  at  time  0.  Some  curves  have  been 
rendered  in  different  line  styles  to  increase  clarity  at  the  intersections. 


do  not  deorbit,  and  constitutes  1.0%  of  the  total  legacy 
damage.  In  contrast,  as  a  near-term  practical  problem,  the 
ASAT  test  can  be  expected  to  cause  the  same  damage  over 
the  next  100  years  as  16.3  rocket  body  insertions  or  41.1 
launches  of  a  spacecraft  that  does  not  deorbit,  and  is 
5.8%  of  the  total  legacy  damage. 

Our  damage  calculations  assume  that  the  deorbit  com¬ 
pliance  rate  9d  and  the  launch  rate  do  not  change  in 
the  future.  If  there  is  a  general  consensus  in  how  these 
parameters  might  change  over  time  (e.g.,  9d  increases  from 
its  current  rate  to  a  rate  where  the  only  noncompliance  is 
due  to  propulsion  system  failures),  then  these  predictions 
could  be  incorporated  into  the  calculations.  With  the 
exception  of  rocket  body  breakups,  the  damage  of  the  var¬ 
ious  activities  varies  by  <12%  for  any  deorbit  compliance 
rate  9d  e  [0, 0.94] .  Complex  interaction  effects  make  dam¬ 
age  calculations  sensitive  to  9d  for  9d  >  0.94. 


3.5.  Allocating  the  legacy  damage 

Finally,  we  consider  the  problem  of  attributing  the  total 
legacy  damage  to  the  various  responsible  parties  (e.g., 
countries).  Let  piy  be  the  proportion  of  objects  of  type 
y  e  Uh  having  party  i  as  their  origin;  we  are  not  aware  of 
any  fragmentation  event  in  which  an  object  from  one 
nation  has  been  fragmented  by  an  object  from  another, 
which  allows  us  to  assume  that  all  fragments  from  an 
object  owned  by  party  i  should  be  attributed  to  it.  As  usual, 
let  y  itself  be  the  number  of  objects  of  type  y.  Finally,  let  coy 
be  the  relative  hazard  of  an  object  of  type  y.  The  propor¬ 
tion  of  the  total  hazard  that  a  party’s  objects  contribute 
is  hi  =  Y^yeuhC°yPiy’  ar|d  the  absolute  hazard  attributed  to 
party  i  is  htL,  where  L  is  the  total  legacy  damage. 

The  proportions  piy  must  be  determined  by  careful  anal¬ 
ysis  of  the  present  catalog,  which  is  beyond  the  scope  of 
this  study.  The  relative  hazard  coy  =  ,  where  the  haz¬ 

ard  wy  is  a  physical  property  of  a  particular  object  type  y 
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relative  to  the  current  orbital  environment,  and  can  be  cal¬ 
culated  as  the  partial  derivative  of  the  damage  nd  with 
respect  to  an  object  type  y: 

~  _dnd  nd(y  +  8y)  -  nd(y) 

dy  dy  ’  1  J 

where  the  argument  of  nd  indicates  what  is  varied.  We  set 
the  perturbation  8y  =  1  (although  the  results  are  insensitive 
to  dy;  0.01  and  1  produce  essentially  the  same  hazards),  so 
that  cby  coincides  with  the  quantity  nd  —  defined  earlier. 
We  consider  an  infinite  horizon  and  find  that  the  relative 
hazard  is  0.239  for  rocket  bodies,  0.739  for  spacecraft  that 
do  not  deorbit,  0.009  for  hazardous  rocket  body  fragments, 
and  0.013  for  hazardous  spacecraft  fragments.  Using  the 
initial  conditions  in  Section  F  and  these  values,  we  find  that 
22.7%,  75.7%,  0.5%,  and  1.1%  of  the  total  legacy  damage  is 
due  to  intact  rocket  bodies,  intact  spacecraft,  hazardous 
rocket  body  fragments,  and  hazardous  spacecraft  frag¬ 
ments,  respectively. 

4.  Discussion 

4.1.  The  model 

Eqs.  (l)-(5)  attempt  to  capture  the  mean  behavior  of  a 
three-dimensional  object-by-object  simulation  such  as 
LEGEND  (Liou  et  al.,  2004)  and  two-dimensional  models 
of  similar  complexity  (e.g.,  SDM  in  (Rossi  et  al.,  1998)). 
Compared  to  most  previous  ODE  models  (Farinella  and 
Cordelli,  1991;  Kessler,  2000;  Rossi  et  al.,  1994)  (a  nota¬ 
ble  exception  is  STAT  in  (Rossi  et  al.,  1998),  which  con¬ 
siders  a  large  number  of  ODEs,  one  for  each  triple  of 
semimajor  axis,  eccentricity,  and  satellite  mass),  ours  is 
distinctive  in  that  we  estimate  parameters  by  taking 
expectations  of  empirical  distributions,  we  introduce  sep¬ 
arate  tracking  of  hazardous  and  benign  fragments,  which 
makes  sense  given  their  different  mean  characteristics,  and 
we  introduce  the  nonuniformity  factor  that  captures  the 
higher  collision  rates  near  the  poles.  Although  our  model 
appears  at  first  glance  to  be  much  simpler  than  the  simu¬ 
lation  model  in  (Liou  et  al.,  2004),  its  complexity  is  in  cal¬ 
culating  parameter  values,  which  are  expectations  over  the 
same  distributions  that  govern  an  object-by-object  simula¬ 
tion.  The  advantages  of  our  approach  are  enumerated 
below. 

A  comparison  to  results  in  Liou  and  Johnson  (2008) 
(which  uses  LEGEND)  shows  that  our  model  generates 
qualitatively  similar  results  (Fig.  1  and  Table  1),  although 
any  number  of  physical  phenomena  and  simple  accounting 
differences  could  explain  the  discrepancies.  These  include: 
differences  in  assumed  characteristics  of  intact  spacecraft 
and  rocket  bodies  (we  assume  all  spacecraft  and  rocket 
bodies  have  the  same  mass  and  take  the  average  area  and 
mass  of  a  satellite  from  Kessler  (2000)  because  Liou  and 
Johnson  (2008)  does  not  provide  the  relevant  quantities); 
different  initial  conditions  (our  knowledge  of  the  initial 


conditions  used  in  Liou  and  Johnson  (2008)  is  quite  lim¬ 
ited);  we  use  mean  values  to  determine  the  effects  of  colli¬ 
sions,  whereas  each  collision  in  Liou  and  Johnson  (2008) 
is  different;  the  cost  of  neglecting  fragments  from  fragmen¬ 
tation  events  occurring  outside  the  SOI  -  though  fragments 
from  above  decay  slowly,  eccentric  fragments  would  imme¬ 
diately  enter  the  SOI;  the  error  made  in  assuming  that 
intact-intact  collisions  occur  according  to  a  chemical  kinet¬ 
ics  model  (in  particular,  intacts  may  have  been  originally 
placed  in  nonconflicting  orbits,  whereas  the  chemical  kinet¬ 
ics  model  assumes  a  particle  is  placed  in  a  trajectory  inde¬ 
pendent  of  the  positions  of  any  of  the  other  particles); 
different  orbital  decay  behavior;  and  rocket  bodies  are 
allowed  to  explode  during  the  first  10  years  of  the  simula¬ 
tion  in  Liou  and  Johnson  (2008). 

There  are  several  advantages  to  using  a  simple  ODE 
model  such  as  Eqs.  (l)-(5)  rather  than  a  complex  simulation 
model:  it  runs  almost  instantaneously,  which  elucidates 
asymptotic  behavior  (Kessler,  2000)  (this  is  particularly 
important  given  our  sustainability  framework,  which 
requires  an  analysis  over  an  infinite  time  horizon)  and  eases 
the  sensitivity  analysis  task  (once  the  parameters  are  calcu¬ 
lated,  one  can  quickly  compute  many  derived  quantities  that 
depend  on  perturbations  in  mean  behavior,  whereas  it  would 
be  necessary  to  run  a  computationally  intensive  object-by- 
object  simulation  repeatedly  to  obtain  mean  behavior,  and 
then  repeat  this  for  each  desired  perturbation);  it  allows  a 
level  of  transparency  and  understanding  that  could  be  valu¬ 
able  if,  e.g.,  multiple  countries  need  to  coalesce  around  a  set 
of  damage  estimates;  and  it  is  amenable  to  rigorous  analyti¬ 
cal  results  involving  asymptotic  behavior.  Our  model  is  gen¬ 
eric  and  can  be  used  at  other  altitudes,  although  the  space 
discretization  may  need  to  be  smaller  at  lower  altitudes  to 
account  for  the  higher  sensitivity  of  escape  rate  with  altitude, 
and  it  may  be  less  accurate  at  altitudes  that  have  consider¬ 
able  congestion  at  shells  just  above  it  (unlike  the  SOI,  see 
Fig.  7  in  Liou  and  Johnson  (2008)). 

Nonetheless,  while  our  qualitative  results  are  likely  to  be 
correct  (e.g.,  the  number  of  fragments  remain  bounded,  the 
sustainable  lifetime  risk  to  an  operational  spacecraft  is  on 
the  order  of  10~3  if  there  is  high  postmission  deorbit  com¬ 
pliance,  and  the  relative  values  of  the  damage  due  to  vari¬ 
ous  activities  maintain  the  rankings  in  Fig.  6),  high-fidelity 
three-dimensional  models  are  required  to  confirm  them. 
Moreover,  the  numerical  output  of  our  model  should  be 
viewed  as  rough  estimates,  and  more  complex  models  are 
required  to  refine  the  numerical  estimates  in  Figs.  2-6. 

4.2.  Results 

Although  our  model  and  analysis  have  some  similarities 
with  previous  research,  and  our  numerical  results  appear  to 
be  similar  to  that  of  previous  research,  we  come  to  different 
conclusions,  both  in  terms  of  analytical  results  and  policy 
implications.  In  terms  of  analysis,  results  in  Kessler 
(2000)  (when  expressed  in  terms  of  our  Eqs.  (6)  and  (7)) 
claim  that  F(t)  grows  without  bound  if  7(f)  >  which  is 
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referred  to  as  the  “runaway  threshold”.  Setting  F(t)  =  0  in 
Eq.  (7)  and  solving  for  F  yields, 

=  (18) 

which  coincides  with  Eq.  (3)  in  Kessler  (2000).  The  discrep¬ 
ancy  between  our  claim  that  F(t)  is  always  bounded  and 
the  result  in  Kessler  (2000)  is  that  the  number  of  intacts  / 
in  the  numerator  of  Eq.  (18)  is  treated  as  a  constant  in 
Kessler  (2000),  whereas  we  capture  via  Eq.  (6)  the  dynamic 
interaction  between  /  and  F  (i.e.,  intacts  decrease  as  frag¬ 
ments  increase). 

More  importantly,  while  our  numerical  results  mimic 
earlier  results  (Liou  and  Johnson,  2005;  Walker  and  Mar¬ 
tin,  2004)  that  stressed  the  importance  of  postmission 
deorbiting,  we  do  not  necessarily  agree  with  the  claim  that 
the  only  way  to  prevent  future  problems  is  to  remove 
existing  large  intacts  from  space  (Liou  and  Johnson, 
2006,  2008).  The  divergence  between  our  views  and  those 
in  Liou  and  Johnson  (2006,  2008)  is  perhaps  due  to  the 
different  performance  metrics  used.  The  root  causes  for 
alarm  in  Liou  and  Johnson  (2006,  2008)  appear  to  be 
the  growth  rate  of  fragments  and  the  small  increase  in 
the  rate  of  catastrophic  collisions  over  the  next  200  years 
(Liou  and  Johnson,  2008,  Fig.  2).  However,  the  great 
majority  of  catastrophic  collisions  in  the  SOI  do  not 
involve  operational  spacecraft,  and  are  hazardous  only 
in  the  sense  that  the  fragments  generated  from  such  a  col¬ 
lision  could  subsequently  damage  or  destroy  operational 
spacecraft.  Therefore,  we  introduced  the  notion  of  the  life¬ 
time  risk  of  an  operational  spacecraft  as  the  primary  per¬ 
formance  metric.  Our  model  predicts  that  the  lifetime  risk 
is  <  5  x  1(T4  over  the  next  two  centuries,  and  always 
stays  <  10-3  if  there  is  very  high  (>98%)  spacecraft  deor¬ 
biting  compliance.  These  risks  appear  to  be  low  relative  to 
the  immense  cost  and  considerable  technological  uncer¬ 
tainty  involved  in  removing  large  objects  from  space,  are 
dwarfed  by  the  «  20%  historical  mission-impacting  (but 
not  necessarily  mission-ending)  failure  rate  of  spacecraft 
(Frost  and  Sullivan,  2004),  and  could  be  overestimated 
if  improved  traffic  management  techniques  lower  future 
collision  risks  (Johnson,  2004).  Hence,  the  need  to  bring 
large  objects  down  from  space  does  not  appear  to  be  as 
clear  cut  as  suggested  in  Liou  and  Johnson  (2006,  2008). 
Nonetheless,  our  model  does  not  incorporate  the  possibil¬ 
ity  of  intentional  catastrophic  collisions  (ASAT  tests, 
space  wars)  that  could  conceivably  occur  in  the  future. 
In  addition,  Fig.  5  considers  only  catastrophic  collisions, 
whereas  noncatastrophic  intact-fragment  collisions  could 
easily  disable  an  operational  spacecraft.  If  the  operational 
lifetime  risk  is  modified  to  include  noncatastrophic  colli¬ 
sions  with  fragments  5=  10cm,  then  the  sustainable  risk 
rises  by  «  50%:  it  increases  from  2.19  x  10-2  to 
3.09  x  10' 2  in  the  base  case,  and  increases  from 
4.91  x  10~4  to  7.94  x  10  4  in  the  full  compliance  case. 
Moreover,  if  fragments  >  1  cm  (rather  than  ^10  cm) 


are  harmful  to  spacecraft  (Johnson,  2004),  then  we  (as 
well  as  other  researchers)  could  be  underestimating  the 
risk. 

In  summary,  in  the  absence  of  the  removal  of  large 
objects  from  space,  the  sustainable  lifetime  risks  in  Figs. 
3-5  do  not  appear  to  be  obviously  above  or  below  a  toler¬ 
able  level.  Even  if  these  risks  are  deemed  acceptable,  it  is 
prudent  to  invest  in  research  and  development  for  space 
remediation  technologies,  which  is  a  topic  of  current  study 
(Proposal  for  forming  an  IAA  study  group,  2000).  How¬ 
ever,  given  the  optimality  of  full  deorbit  compliance  from 
a  societal,  sustainable  perspective,  and  the  sensitivity  of 
sustainable  lifetime  risk  to  postmission  deorbit  compliance, 
the  primary  focus  for  policymakers  should  be  on  increasing 
compliance,  which  leads  us  to  a  discussion  of  economic 
instruments  that  could  be  used  to  address  this  issue. 

4.3.  Setting  fees 

Creating  a  set  of  viable  economic  instruments  will  be 
extremely  difficult,  and  requires  sustained  political  negotia¬ 
tions  by  a  large  number  of  actors;  here,  we  simply  discuss 
some  of  the  key  challenges  that  would  be  faced  in  trans¬ 
forming  the  damages  in  Fig.  6  into  specific  fees  that  could 
be  imposed  on  past  and  future  debris-generating  activities. 
Such  fees  are  consistent  with  the  “polluter  pays”  principle 
proposed  to  address  global  climate  change  (Claussen  et  al., 
1998). 

At  a  strategic  level,  the  first  challenge  is  to  decide  what 
types  of  activities  to  tax,  and  for  what  reasons.  One  reason 
to  charge  a  fee  is  to  deter  future  debris-generating  activity 
(i.e.,  ideally  no  fees  are  collected  for  certain  activities), 
while  another  reason  could  be  to  generate  funds  for,  e.g., 
compensating  owners  of  operational  spacecraft  that  are 
victims  of  future  catastrophic  collisions  (which  are  largely 
random  events),  subsidizing  postmission  deorbiting,  or  per¬ 
forming  research  and  development  into  various  mitigation 
technologies  (e.g.,  methods  that  remove  large  objects  from 
space,  electric  propulsion  systems  that  could  reduce  the 
cost  of  deorbiting  (Ryden  et  al.,  1997),  improved  space 
traffic  management  techniques  (Johnson,  2004),  and  space¬ 
craft  shielding). 

We  begin  with  the  most  benign  activity,  which  is  the 
launching  of  a  spacecraft  that  deorbits  after  its  mission. 
Because  this  spacecraft  is  doing  its  best  to  preserve  space, 
one  might  argue  that  it  should  only  be  charged  a  fee  that 
compensates  for  its  future  expected  damage.  In  this  case, 
the  damage  that  occurs  at  time  t  should  be  multiplied  by 
e-rt  inside  the  integral  in  Eq.  (16),  where  r  is  the  discount 
rate,  and  the  total  discounted  damage  should  be  multiplied 
by  the  replacement  cost  of  a  satellite.  As  an  illustrative 
example,  if  we  set  r  =  0.05  (Perman  et  al.,  2003)  and 
Cs  =  $0.5B,  then  the  launch  fee  is  only  $980  (for  implemen¬ 
tation,  the  values  of  r  and  Cs  need  to  be  estimated  with 
considerable  care  and  achieve  broad  consensus). 

Note  that  a  noncompliance  fee  for  spacecraft  that  do 
not  deorbit  can  be  viewed  as  using  differential  launch  fees 
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that  depend  on  whether  or  not  a  spacecraft  deorbits.  The 
primary  motivation  of  a  noncompliance  fee  is  to  preserve 
space  (Fig.  5).  To  produce  a  deterrent  effect,  the  fee  to 
launch  a  spacecraft  that  does  not  deorbit  needs  to  be  at 
least  the  deorbit  cost  plus  the  fee  (if  any)  to  launch  a  space¬ 
craft  that  deorbits.  The  voluntary  guideline,  which  implies 
that  it  is  better  to  incur  a  deorbit  cost  now  rather  than 
destroy  0.065  spacecraft  in  the  future,  according  to  our  cal¬ 
culations,  has  not  achieved  full  compliance.  It  is  important 
to  understand  the  reasons  for  the  current  deorbit  noncom¬ 
pliance,  which  could  shed  light  on  how  to  improve 
compliance.  If  players  are  acting  rationally  and  in  their 
own  self-interest,  then  possible  reasons  for  not  deorbiting 
are  (i)  they  are  behaving  myopically  and  are  discounting 
collisions  that  happen  far  in  the  future,  perhaps  having 
faith  in  a  future  technological  solution  to  the  problem; 
(ii)  they  do  not  plan  many  future  launches  in  congested 
regions  of  space  (such  as  the  SOI),  and  so  will  not  receive 
the  benefits  (i.e.,  avoided  future  collisions)  from  deorbiting; 
and  (iii)  they  have  an  inventory  of  old  spacecraft  that  are 
expensive  to  retrofit  for  deorbiting.  To  the  extent  that  dif¬ 
ferent  players  have  different  (perceived)  deorbit  costs  and 
benefits,  the  noncompliance  fee  (or  deorbit  subsidy)  needs 
to  be  sufficiently  high  to  deter  all  players.  If  the  fee  is  to 
be  based  on  damages  and  if  a  single  player  (i.e.,  space 
agency)  can  substantially  alter  6d  by  his  postmission  deor¬ 
bit  policy,  then  it  may  be  more  appropriate  to  assess  the 
damage  of  his  aggregate  policy  decision. 

Another  challenge  for  a  noncompliance  fee  stems  from 
the  moral  hazard  issue:  the  failure  rate  of  propulsion  sys¬ 
tems  has  been  estimated  to  be  3.9%  (de  Week  et  al., 
2003,  Fig.  8)  and  hence  even  if  a  spacecraft  owner  plans 
to  deorbit  after  the  mission  is  complete,  he  may  fail  to  do 
so,  and  the  international  space  community  may  have  no 
way  of  verifying  whether  the  failure  to  deorbit  was  deliber¬ 
ate  or  unintentional. 

The  fee  structure  would  have  to  address  this  issue.  Much 
of  the  above  discussion  about  a  deorbit  noncompliance  fee 
also  holds  for  failing  to  deorbit  rocket  bodies  or  failing  to 
follow  passivation  techniques,  although  compliance  for 
these  activities  in  the  SOI  appears  to  be  much  closer  to 
100%  than  spacecraft  deorbiting. 

There  is  also  a  legacy  issue,  which  typically  arises  when 
discussing  launch  costs  (e.g.,  Prasad,  2005),  in  that  most  of 
the  satelhtes  in  LEO  were  launched  by  a  small  group  of  coun¬ 
tries.  Hence,  there  could  also  be  a  one-time  legacy  fee  due  to 
all  existing  space  objects  currently  in  the  SOI.  The  use  of  leg¬ 
acy  fees  is  a  particularly  thorny  political  issue  related  to,  e.g., 
the  Cold  War  and  the  original  research  and  development 
costs  incurred  by  the  U.S.  and  U.S.S.R.  that  subsequent 
countries  did  not  bear.  Sustainability  is  a  modem  goal  and 
during  the  Cold  War  (when  much  of  the  existing  debris 
was  generated),  the  U.S.  and  U.S.S.R.  were  not  concerned 
with  preserving  space  for  the  future.  It  is  also  not  clear  to 
what  extent  Russia  feels  responsible  for  past  U.S.S.R. 
behavior.  In  any  case,  these  activities  are  in  the  past  and 
not  deterable.  If  one  wanted  to  charge  a  one-time  legacy  cost 


to  cover  future  expected  damages,  then  (using  r  =  0.05  and 
Cs  =  $0.5B)  the  legacy  fee  is  S5.55M,  which  is  negligible 
compared  to  the  annual  budgets  of  the  U.S.  and  Russian 
space  programs.  To  allocate  this  fee  across  various  players, 
the  relative  hazards  (cby)  should  also  be  computed  with 
r  =  0.05,  which  yields  0.618  for  rocket  bodies,  0.236  for 
spacecraft  that  do  not  deorbit,  0.077  for  hazardous  rocket 
body  fragments,  and  0.068  for  hazardous  spacecraft  frag¬ 
ments;  combining  these  numbers  with  the  initial  conditions, 
we  find  that  the  total  legacy  damage  due  to  these  four  satel¬ 
lite  types  is  63.0%,  26.0%,  4.6%,  and  6.4%,  respectively. 

ASAT  tests  are  another  politically-charged  issue.  The 
FengYun  1C  ASAT  test  may  have  been  partially  motivated 
by  former  U.S.  and  Soviet  ASAT  tests  and  by  the  U.S. 
December  2001  withdrawal  from  the  Anti-Ballistic  Missile 
treaty.  As  with  a  deorbit  noncompliance  fee,  the  primary 
goal  of  an  ASAT  test  fee  would  be  as  a  deterrent.  In  con¬ 
trast  to  the  case  where  there  is  a  failure  to  deorbit,  it  is  dif¬ 
ficult  to  determine  the  magnitude  of  a  financial  deterrent 
because  it  is  very  difficult  to  put  a  financial  cost  on  the 
value  (to  the  tester)  of  an  ASAT  test.  The  ASAT  test  fee 
might  be  set  to  significantly  exceed  its  damage  so  as  to  dis¬ 
courage  the  militarization  of  space,  essentially  charging  it 
not  only  for  its  test  but  for  the  escalation  effect  it  may  have 
on  other  players.  In  any  case,  a  test  fee  would  be  very  dif¬ 
ficult  to  implement  from  a  political  standpoint. 


5.  Conclusions 

It  appears  that  if  full  compliance  of  the  25-year  space¬ 
craft  deorbiting  guidelines  can  be  achieved  within  the  next 
few  decades  and  no  ASATs  are  used  or  tested  in  the  future, 
then  the  lifetime  risk  from  space  debris  in  the  SOI  may  be 
sustainable  at  a  tolerable  («  10-3)  level.  Hence,  the  focus 
of  policy  should  be  on  achieving  full  deorbit  compliance 
(including  improving  the  reliability  of  propulsion  systems) 
and  fostering  a  culture  that  makes  ASAT  tests  (and  use) 
taboo.  Indeed,  at  the  current  launch  rate  of  3/year  and  a 
deorbit  cost  of  S0.5M,  the  total  undiscounted  deorbit  cost 
after  1000  years  is  S1.5B,  which  is  the  cost  of  launching  a 
single  defense  spacecraft.  It  seems  improbable  that  a  future 
technology  will  be  able  to  clean  up  space  for  less  cost  than 
launching  a  single  defense  spacecraft.  Our  analysis  also 
provides  a  framework  for  setting  fees  for  past  and  future 
space  activities,  but  many  challenges  remain  in  developing 
a  viable  set  of  economic  instruments.  Finally,  we  reiterate 
that  our  analytical  approach  does  not  produce  results  that 
are  as  accurate  as  a  high-fidelity  three-dimensional  model 
(e.g.,  Liou  et  al.,  2004).  Hence,  our  results  need  to  be  con¬ 
firmed  by  a  more  complex  model  before  they  can  inform 
major  space  policy  decisions. 
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Appendix  A 

In  this  appendix,  we  obtain  values  for  the  parameters  in 
Eqs.  (1)— (5);  a  more  detailed  description  of  our  procedure  - 
sufficient  for  reproducibility  of  our  results  -  appears  in  Sec¬ 
tion  B  of  Bradley  and  Wein  (2008).  The  values  for  the  col¬ 
lision  rate  parameters  flay,  the  fragment  generation 
parameters  <5™,  and  the  decay  rates  are  found  by  physi¬ 
cal  modeling  (Alberty  and  Silbey,  1997;  King-Hele,  1964) 
and  calculations  using  the  empirical  probability  distribu¬ 
tions  and  formulas  describing  satellite  fragmentation  and 
fragment  characteristics  developed  in  Johnson  et  al. 
(2001),  Kessler  (2000),  and  Rossi  et  al.  (1994).  The  compu¬ 
tational  framework  for  incorporating  the  data  in  Johnson 
et  al.  (2001)  is  described  in  Section  B.  We  estimate  the  col¬ 
lision  rate  parameters  in  Section  C,  the  fragment  genera¬ 
tion  parameters  in  Section  D,  the  decay  rates  in  Section 
E,  and  the  initial  conditions  in  Section  F. 

Appendix  B.  Computational  framework 

Calculations  for  the  collision  rates,  fragment  generation 
rates  and  decay  rates  in  Sections  C,  D  and  E  require  inte¬ 
grating  over  five  probability  distributions,  which  are  listed 
in  the  next  paragraph,  to  obtain  expected  values.  More¬ 
over,  for  the  collision  rates  and  fragment  generation  rates, 
we  are  interested  in  a  weighted  mean  where  we  associate 
two  weights  with  a  fragment.  The  first  weight  function, 
Wi,  accounts  for  the  proportion  of  time  the  fragment 
spends  in  the  SOI  as  well  as  considerations  related  to  eccen¬ 
tricity  and  deorbit  time.  The  second  weight  function,  w2, 
quantifies  the  hazard  a  fragment  poses  to  an  intact,  and 
essentially  defines  the  fragments  of  interest:  those  that  are 
hazardous  to  intacts. 

Of  the  five  random  variables  that  we  consider,  the 
empirical  probability  distributions  in  Johnson  et  al. 
(2001)  govern  three  of  them: 

Xc  =  log10Zc, 

where  Z,cis  the  characteristic  length,  (B.l) 

X  =  l°gio  iA/M), 

where  A/M  is  the  area-to-mass  ratio,  (B.2) 

v  =  log10|Ai>|, 

where  Av  is  the  change  in  velocity.  (B.3) 

These  distributions  depend  on  whether  the  fragment  ejec¬ 
tion  event  is  a  collision  or  an  explosion  and  whether  the 
original  intact  is  a  rocket  body  or  a  spacecraft.  We  are 


interested  only  in  collisions.  We  introduce  two  additional 
random  variables: 

u,  the  altitude  of  the  satellite  in  circular 

orbit  generating  the  fragments,  (B.4) 

z,  the  direction  of  Av  on  the  sphere.  (B.5) 

We  use  the  same  symbol  for  both  a  random  variable  and  a 
value  drawn  from  its  distribution  if  the  distinction  is  clear 
from  the  context.  The  characteristic  length  Lc  in  (B.l)  has  a 
power  law  distribution,  x  in  (B.2)  is  conditioned  on  Xc  and 
is  a  mixture  of  two  Gaussians,  v  in  (B.3)  is  conditioned  on  x 
and  has  a  normal  distribution,  u  is  independent  and  uni¬ 
formly  distributed  on  the  interval  [900, 1000]  km,  and  z  is 
independent  and  uniformly  distributed  on  the  sphere.  The 
random  vector  Q  has  as  elements  the  five  random  variables 
in  Eqs.  (B.1)-(B.5). 

Later  in  this  section,  we  also  use  Eq.  (9)  from  Johnson 
et  al.  (2001),  which  states  that  the  average  area  of  a  frag¬ 
ment  having  characteristic  length  Lc  is 

A  =  aL2+£,  (B.6) 

where  a.  =  0.556945  and  e  =  0.0047077.  Rather  than  view 
A  as  a  conditional  expectation  given  a  value  of  Lc,  we  treat 
Eq.  (B.6)  as  a  deterministic  mapping  from  the  realization  of 
a  random  characteristic  length  Lc  to  a  realization  of  a  ran¬ 
dom  area  A.  Similarly,  a  fragment  that  is  randomly  as¬ 
signed  the  value  x  from  Eq.  (B.2)  has  a  random  mass  m 
given  by 


where  A  is  the  area  in  Eq.  (B.6).  That  is,  a  fragment’s  mass 
is  determined  by  its  characteristic  length  and  its  A/ M  ratio 
via  Eqs.  (B.l),  (B.2),  (B.6),  and  (B.7). 

The  remainder  of  this  section  is  organized  as  follows:  the 
proportion  of  time  a  fragment  spends  in  the  SOI  is  calcu¬ 
lated  in  Section  B.l  and  the  probability  a  collision  is  cata¬ 
strophic  is  given  in  Section  B.2;  these  two  calculations  form 
the  basis  of  our  weights  wi  and  w2,  respectively.  Finally,  the 
expectation  of  a  generic  function  of  the  random  vector  12  is 
given  in  Section  B.3. 

B.l.  The  fraction  of  time  spent  in  the  SOI 

We  use  classical  orbital  mechanics  equations  (in  partic¬ 
ular,  Eqs.  (3.3),  (3.10),  (3.17),  and  (3.20)  in  King-Hele 
(1964))  and  integrate  the  equations  of  motion  to  deter¬ 
mine  the  fraction  of  time  each  fragment  (whose  properties 
have  been  sampled  via  Eqs.  (B.1)-(B.5))  spends  in  the  SOI 
(see  Section  B.l.l  of  Bradley  and  Wein  (2008)  for  details); 
each  intact  in  our  model  is  assumed  to  be  in  a  circular 
orbit. 

The  weight  w i  (introduced  in  Section  B)  is  defined  to  be 
the  calculated  fraction  of  time  spent  in  the  SOI,  which  var¬ 
ies  according  to  the  fragment’s  properties  in  Eqs.  (B.l)- 
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(B.5).  We  also  use  this  calculation  of  wi  to  help  determine 
the  launch  rate  and  the  initial  number  of  intacts. 

As  an  aside,  we  conclude  this  subsection  by  performing 
a  calculation  (similar  to  one  in  Rossi  (2002))  that  is 
required  in  Section  3.3:  to  determine  the  change  in  velocity 
necessary  to  deorbit  an  intact  from  the  SOI.  We  assume  the 
intact  has  an  area-to-mass  ratio  of  0.01  nr/kg  (Kaula, 
1983)  and  does  not  deploy  an  end-of-life  sail-like  device 
to  increase  this  ratio.  The  intact  is  initially  in  a  circular 
orbit  at  950  km  and  its  initial  velocity  is  7.375  km/s  (Eq. 
(3.17)  in  King-Hele  (1964)).  Fig.  6-2  of  NASA  Safety  Stan¬ 
dard  1740.14  (1995)  gives  orbits  as  (perigee,  apogee)  pairs 
for  various  area-to-mass  ratios  that  decay  naturally  in 
25  years.  We  shall  be  interested  in  two  pairs:  (500,950) 
and  (510,900). 

By  considering  an  apogee  of  950  km  and  a  perigee  of 
500  km,  we  find  that  the  method  that  simply  lowers  the 
perigee  to  500  km  requires  a  speed  at  apogee  of 
7.257  km/s.  Hence,  the  required  change  in  velocity  is 
Av  =  118  m/s.  A  second  method  first  lowers  the  perigee 
to  510  km  and  then  lowers  the  apogee  to  900  km,  and 
requires  Av  =  128  m/s.  The  first  method  lets  the  spacecraft 
linger  in  the  SOI  for  a  little  while;  the  second  method  shows 
that  only  10.8%  more  fuel  lets  the  spacecraft  leave  the  SOI 
essentially  immediately.  Our  model  -  by  assuming  all  deor¬ 
bit-capable  spacecraft  are  operational  and  without  provid¬ 
ing  a  nonoperational  deorbit  state  -  implicitly  assumes  the 
second  method  is  used,  although  our  estimate  in  the  main 
text  of  the  deorbit  cost  as  a  function  of  Av  is  too  crude 
to  distinguish  between  the  two  methods. 


B.2.  The  probability  a  collision  is  catastrophic 

In  the  model  in  Johnson  et  al.  (2001),  whether  a  collision 
is  catastrophic  depends  on  the  relative  velocity  and  masses 
of  the  objects.  Because  we  assume  that  rocket  bodies  and 
spacecraft  have  the  same  mass,  a  particular  projectile  is 
hazardous  to  these  two  satellites  with  equal  probability. 

Fig.  9  of  Johnson  et  al.  (2001)  presents,  for  altitudes  of 
200,  500,  1000,  and  1500  km,  a  piecewise-linear  probability 
density  function  fVc(vc)  of  the  relative  velocity  of  two 
objects  given  that  they  collide  (i.e.,  the  relative  collision 
velocity).  We  use  the  distribution  for  the  altitude 
1000  km,  and  hence  the  result  is  slightly  conservative  rela¬ 
tive  to  what  we  would  obtain  using  a  distribution  for  the 
altitude  950  km,  which  is  the  center  of  our  SOI.  This  distri¬ 
bution  implies  that  the  mean  relative  collision  velocity  is 

vc  =E[vc]  «  12.1  km/s,  (B.8) 

which  we  shall  need  later. 

Fet  two  colliding  objects  have  masses  mi  and  m2  >  mt 
(in  kg).  According  to  Johnson  et  al.  (2001),  the  collision 
is  catastrophic  if  the  kinetic  energy  of  the  projectile  in 
Joules  is  40  times  greater  than  the  mass  of  the  target  in 
grams.  Hence,  there  is  a  critical  velocity  depending  on  m, 
and  m2,  which  we  denote  by  rcrit(mi .  w2),  such  that 


vc  ^  rCrit(wi , mi)  results  in  a  catastrophic  collision.  Then 
the  probability  the  collision  is  catastrophic  is 

pc{mum2)  =  J^  ^  ^fVc(v)dv.  (B.9) 

A  fragment’s  w2  weight  is  the  probability  that  an  intact- 
fragment  collision  involving  the  given  fragment  will  be  cat¬ 
astrophic.  Kessler  (2000)  assumes  that  any  combination  of 
two  intacts  have  mass  1600  kg,  and  we  assume  that  every 
intact  has  a  mass  of  M  =  800  kg.  In  our  simulation  model, 
each  fragment  is  assigned  a  random  mass  via  Eq.  (B.7).  A 
hazardous  fragment  ( Fhz )  with  mass  m  has  weight 
w2  =  pc(m,M),  and  a  benign  fragment  (Fb)  with  mass  m 
has  weight  w2  =  1  -  pc(m,  M),  where  pfm,  M)  is  derived 
from  Eq.  (B.9).  In  Section  C,  we  also  use  Eq.  (B.9)  to  cal¬ 
culate  the  fractions  of  intact-intact  and  fragment-fragment 
collisions  that  are  catastrophic. 


B.3.  Computing  expectations 


Because  wi  depends  on  the  fragment’s  orbit,  it  is  a  func¬ 
tion  of  v,  u,  and  z  in  Eqs.  (B.3)-(B.5).  Moreover,  v  is  con¬ 
ditioned  on  x,  and  y  is  conditioned  on  lc.  Hence,  any 
expectation  involving  the  weighting  function  wi  must  inte¬ 
grate  over  all  five  probability  density  functions  (PDFs). 
The  weight  w2  is  a  function  of  the  mass  m  and  the  hazard 
type  KG{h,b};  we  write  w2(-;k)  to  indicate  the  latter 
dependence.  The  mass  m  is  a  function  of  the  area  A  and 
X  via  Eq.  (B.7),  and  A  is  a  function  of  Xc  by  Eq.  (B.l).  Addi¬ 
tionally,  the  PDFs  of  some  of  the  random  variables  are 
parameterized  by  the  source  of  the  fragments  (i.e.,  rocket 
body  or  spacecraft). 

Let  /n(o);r,  k)  be  the  PDF  of  the  random  vector  Q 
parameterized  by  the  source  r  e  U1  and  the  hazard  type 
k  e  {h,b},  and  define  the  weights 


w(co-,  k)  = 


wi(co)w2((q;k) 
Jwi(co)w2(co-,k){  d"' 


(B10) 


Consider  first  a  generic  function  g(oj).  The  weighted 
expected  value  of  g  is 

Ew[g;  T,  k]  =  J  g{(Q)fa{w,X,K)w{w,K)&(Q.  (B.ll) 


Similarly,  consider  a  generic  function  h{ao\,  m2),  a  function 
of  two  vectors  independently  drawn  from  the  distribution 
fa.  The  weighted  expected  value  of  li  is 

Eww[h;xi,KUx2,K2\  =  J  h((oU(o2)fc2(coi;TUKi) fa 

x  ( (o2;  x2,  K2)w(a>i-,  Ki)w(a>2,  tc2) dwi  dco2 .  (B.12) 


In  Sections  C,  D  and  E,  we  compute  the  collision  rates,  the 
fragment  generation  rates,  and  the  decay  rates  by  substitut¬ 
ing  specific  functions  for  g(co)  and  h(a>i,co2)  using  Eqs. 
(B.ll)  and  (B.12). 

We  approximate  the  expected  values  in  Eqs.  (B.ll)  and 
(B.12)  by  sampling  the  random  vector  Q  repeatedly  and 
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computing  the  sample  weighted  mean  of  the  function  at 
each  sample  point;  details  can  be  found  in  Section  B.1.3 
of  Bradley  and  Wein  (2008). 

Appendix  C.  Collision  rate  parameters 

We  derive  an  analytical  expression  for  the  collision  rate 
parameter  (1  in  Section  C.l,  estimate  the  parameter  values 
in  this  analytical  expression  in  Section  C.2,  and  incorporate 
nonuniform  spatial  density  in  the  SOI  in  Section  C.3. 

C.l.  Collisions  between  spherical  objects 

We  assume  satellite  collisions  can  be  modeled  by  the 
ideal-gas  model  in  chemical  kinetics  (Alberty  and  Silbey, 
1997).  Consider  two  spheres  contained  in  a  large  volume 
V having  random  radii  n  and  r2,  and  random  relative  speed 
On  ~  fvjv).  They  collide  if  their  centers  are  within  a  dis¬ 
tance  d]2  =rj+  r2,  referred  to  as  the  collision  diameter, 
of  each  other.  As  in  Johnson  et  al.  (2001)  and  Kessler 
(2000),  we  use  the  collision  cross  section  o12  rather  than 
the  collision  diameter  dn,  where  =  nd2n.  Let 
f  12  =  / vu fVi2(v)  dr  be  the  collision  probability  velocity.  It 
follows  (see  Section  B.2. 1  of  Bradley  and  Wein  (2008)  for 
details)  from  Eqs.  (l)-(5)  that  the  parameter  (l  corresponds 
to  ’fE[on]. 

Finally,  while  the  w2  weights  in  Section  B.2  quantify  the 
fraction  of  intact-fragment  collisions  that  are  catastrophic, 
we  have  yet  to  account  for  the  fraction  of  intact-intact  and 
fragment-fragment  collisions  that  are  catastrophic.  We  do 
so  using  the  variable  where  a  proportion  \  of  colli¬ 
sions  are  disregarded  because  they  are  not  catastrophic;  as 
shown  in  Section  C.2,  t j/  is  derived  for  various  intact-intact 
and  fragment-fragment  collisions  using  the  probabilities 
pc(mum2)  in  Eq.  (B.9).  Hence,  we  define  the  collision  rate 
parameter  by 

P  =  V-fE[^n],  (C.l) 

where  the  tilde  denotes  that  the  spatial  density  of  satellites 
is  assumed  to  be  uniform  in  the  SOI.  In  Section  C.3,  we 
multiply  /?  by  a  nonuniformity  factor  to  obtain  the  collision 
rate  parameter  /l. 

C.2.  Estimating  parameters  in  Eq.  (C.l) 

In  this  subsection,  we  estimate  the  quantities  V,  v12,  <*n, 
and  i ji  in  Eq.  (C.l).  The  Earth’s  radius  is  6378  km,  and 
hence  the  volume  of  our  compartment  is  V  =  |7t((6378+ 
1000)3  -  (6378  +  900)3)  =  6.748  x  1010  km3.  We  take  the 
collision  probability  velocity  to  be  vu  =  8.8  km/s  (Johnson 
et  al.,  2001,  Fig.  8). 

For  a. ,ye{t/},  we  have  ^/ojj  =  sfa)  +  y/o),  and  it 
remains  to  determine  the  cross  sections  os,oR,  and  o>* .  In 
his  calculations,  Kessler  (Kessler,  2000)  assumes  that  the 
collision  cross  section  of  a  fragment  is  zero;  hence  his  cross 
section  values  are  for  intacts  only.  We  do  not  know  the 


individual  cross  sections  of  spacecraft  or  rocket  bodies 
and  so  must  use  average  values  to  obtain  as  and  aR.  Kessler 
obtains  average  cross  sections  for  the  region  [920, 1020]  km 
of  <7/=  14.1m2  for  intact-fragment  collisions  and 
g i  =  53.8m2  for  intact-intact  collisions  (Kessler,  2000). 
The  ratio  cv/cy  =  3.82  <  4  (by  less  than  5  %)  because  c, 
is  computed  by  averaging  over  all  intact-intact  combina¬ 
tions  in  his  database.  Our  calculations  assume  the  ratio  is 
exactly  4,  effectively  using  only  the  14.1  m2  value.  Kessler 
observes  that  spacecraft  have  an  average  cross  section 
0.35  times  the  average,  or  as  =  0.35  x  14.1  =  4.94  m2, 
and  rocket  bodies  have  an  average  cross  section 
oR  =  1.65  x  14.1  =  23.3  m2  (Kessler,  2000). 

For  fragments,  we  have  a  probability  distribution  over 
the  area  A  from  Eq.  (B.6),  which  corresponds  to  the  cross 
section  under  the  assumption  that  fragments  are  spherical. 
Under  this  assumption,  the  collision  cross  section  for  an 
intact  and  fragment  is 

=  ^\/^a  +  ^/aJ  ,  (C.2) 

and  the  collision  cross  section  for  two  fragments  is 

aCx\F^  =  (V*  “I"  •  (C.3) 

The  final  parameters  to  estimate  are  the  catastrophic 
probabilities  i/f.  The  change  over  time  of  the  proportion 
of  fragments  that  are  hazardous  to  intacts  is  accounted 
for  by  the  separate  hazardous  and  benign  fragment  types 
and  corresponding  equations  in  Eqs.  (l)-(5)  and  the  w2 
weights  in  Section  B.2.  Therefore,  all  intact-fragment  colli¬ 
sions  are  counted  and  we  set  \j/xF*  =  1  for  a  e  U1. 

All  intacts  have  the  same  mass,  M  =  800  kg.  Therefore, 
each  poses  the  same  threat  to  the  others,  and  Eq.  (B.9) 
implies  that 

t/y  =  t/^  =  pc(M,M)  =  0.9996  for  a,  y  €  f/7;  (C.4) 

i.e.,  nearly  all  intact-intact  collisions  are  catastrophic. 

Fragments  pose  different  hazards  to  each  other.  By  Eq. 
(B.9),  the  probability  that  a  collision  between  a  fragment  of 
type  FKz'x  with  mass  mx  and  one  of  type  F((  with  mass  m2  is 
catastrophic  is 

i/V*i  =  pc(mi ,  m2) .  (C.5) 

For  all  types  of  fragment-fragment  collisions,  1/7  >  0.96 
when  averaged  over  Q. 

Though  Johnson  et  al.  (2001)  defines  a  noncatastrophic 
collision  as  one  in  which  only  the  less  massive  object  breaks 
up,  we  assume  that  in  intact-intact  and  fragment-fragment 
collisions,  a  noncatastrophic  collision  produces  no  new 
fragments.  Because  essentially  all  intact-intact  and  almost 
all  fragment-fragment  collisions  are  catastrophic,  there  is 
little  error  in  making  this  assumption,  and  the  error  is  con¬ 
servative.  As  in  Johnson  et  al.  (2001),  we  assume  that  only 
the  less  massive  object  breaks  up  in  noncatastrophic  intact- 
fragment  collisions. 
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The  value  of  fi  is  obtained  by  taking  the  expectation  in 
Eq.  (C.l)  and  substituting  Eqs.  (C.2)-(C.5)  for  the  generic 
functions  in  Eqs.  (B.ll)  and  (B.12).  For  intact-intact  colli¬ 
sions,  Eq.  (C.l)  contains  no  random  variables  and  we  have 

Pay  =  f°r  (C.6) 

where  t/q  is  given  in  Eq.  (C.4).  For  intact-fragment  colli¬ 
sions,  we  have 

^=y£w[^;r,*c],  (C.7) 

and  for  fragment-fragment  collisions,  we  have 

[i/h^i f* ffF*; ;  Ti,Ki,t2,  k2 j .  (C.8) 

The  numerical  values  from  Eqs.  (C.6)-(C.8)  for  both  ver¬ 
sions  of  the  model  (i.e.,  in  the  absence  and  presence  of  frag¬ 
ment-fragment  collisions)  appear  in  Tables  C.l  and  C.2. 

C.3.  Nonuniform  spatial  density 

Satellites  are  distributed  nonuniformly  throughout  the 
SOI  for  two  reasons:  a  satellite  having  an  approximately 


circular  orbit  tends  to  spend  a  greater  proportion  of  time 
at  latitudes  near  its  inclination  than  at  other  latitudes, 
and  the  volume  of  the  SOI  as  a  function  of  latitude 
decreases  near  the  poles.  We  incorporate  these  phenomena 
into  our  model  by  considering  a  nonuniform  version  of  the 
ideal  gas  model. 

Consider  populations  of  two  types  of  spherical  objects 
mixed  homogeneously  in  a  region  of  space  having  volume 
V.  Each  population  has  number  density  pt  =  f,i—  1,2. 
The  ideal  gas  model  gives  the  number  of  collisions  c  in  a 
unit  of  time  as 

c  =  bnVpxp2  =  bnn^  (C.9) 

for  some  constant  b]2. 

Suppose  the  objects  prefer  some  parts  of  space  over  oth¬ 
ers.  Break  the  region  into  voxels  indexed  by  j,  each  having 
volume  V j.  Let  the  objects  of  type  i  have  number  density 
HiPij  in  voxel  j.  Eq.  (C.9)  holds  in  each  voxel  j;  summing 
over  the  voxels  yields  c  =  b\2n\n2Y^jV jPyPij  =  fubn^jf, 
where  we  define 

fn  =  VY^VjPyPy  (C.10) 


Table  C.l 

Parameter  values  for  Eqs.  (1)— (5),  where  time  is  in  years.  In  each  of  the  five  sections  of  the  table,  the  parameter  is  given  in  the  far  left  column,  the  rows  are 
the  first  subscript,  the  columns  are  the  second  subscript,  and  the  table  entries  are  the  parameter  values;  e.g.,  the  second  row,  fifth  column  of  the  second 
section  yields  =  -1.00. 
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Table  C.2 

Parameter  values  for  the  fragment-fragment  version  of  Eqs.  (l)-{5).  See  the  legend  to  Table  C.l  for  an  explanation. 
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fr 

fr 
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-0.04 

fr 

0 

17.87 

1.35 

FS 

41.53 

21.27 

Fs 

2.94 

to  be  the  nonuniformity  factor.  If  the  density  is  uniform, 
then  riiPy  =  f  and  hence  fn=  1 . 

We  assume  that  the  density  ptj  is  constant  in  time,  so 
that  c(t)  =  fnb\2n\{t)n2{t) / V  for  all  t  >  0.  Hence,  we  can 
accommodate  spatial  nonuniformity  in  our  ODE  model 
by  simply  multiplying  the  rate  /?  by  the  nonuniformity  fac¬ 
tor.  More  specifically,  in  the  remainder  of  this  section  we 
calculate  fay  for  ac,  y  e  {R.S,  F},  which  allows  us  to  use 
Eqs.  (C.6)-(C.8)  to  define  the  collision  rate  parameters 

P*y  =  .UA-,  for  a,  y  €  U\  (C.  1 1 ) 

/U«  =  Uhn  for  *  6  U',  (C.12) 

(C.13) 

We  discretize  SOI  into  90  voxels,  where  voxel 
j  6  {1.2, ... ,  90}  spans  the  latitudes  (J  -  1)°  to  j°,  so  that 
the  volume  of  voxel  j  is 

/*2  I  /so  r^e+ 1000  r2n 

Vj=  /  /  r2  sin  (j)  dGdrdcj).  (C.14) 

“2^~^180  d  Re+ 900  Jo 

For  each  fragment,  rocket  body,  and  spacecraft  in  the  cur¬ 
rent  publicly  available  catalog  (Space  Track,  2008),  we  use 
the  method  described  in  Section  B.l  to  calculate  the  pro¬ 
portion  of  time  the  satellite  spends  in  the  SOI  in  each  lati- 


Table  C.3 

The  nonuniformity  parameters  fay  derived  in  Section  C.3.  The  rows 
correspond  to  the  first  subscript  and  the  columns  correspond  to  the  second 
subscript;  e.g„  fm  =  1.39. 


/ 

F 

R 

S 

F 

1.38 

1.39 

1.33 

R 

1.72 

1.55 

S 

1.44 

tude  bin  (see  Section  B.2.3  of  Bradley  and  Wein  (2008)  for 
details),  which  allows  us  to  determine  the  densities  pir  The 
nonuniformity  factors  fxy  are  calculated  via  Eq.  (C.10)  and 
appear  in  Table  C.3.  An  analysis  of  the  catalog  shows  that 
there  are  a  large  number  of  SL-8  rocket  bodies  having  incli¬ 
nations  between  82°  and  83°;  their  presence  is  the  primary 
reason  why  fRy  is  larger  than  the  other  nonuniformity 
factors. 

Appendix  D.  Debris  generation  parameters 

In  this  section,  we  derive  an  expression  for  the  debris 
generation  parameter.  Suppose  that  in  addition  to  a  radius 
R,  a  sphere  has  a  mass  M,  and  that  the  random  variables 
(. R,M )  =  Q  ~  /~(w).  Colliding  spheres  break  up  to  pro- 
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duce  fragments.  A  collision  between  spheres  1  and  2  gener¬ 
ates  N(mi,m2 )  fragments.  Then  the  expected  number  of 
fragments  generated  per  unit  time  is 

j Ar(mi,/n2)(r1+r2)2/~(d;i)/~(d)2) dxo1da2.  (D.l) 

By  Eqs.  (1)— (5),  we  see  that  the  debris  generation  parame¬ 
ter  <5  corresponds  to 

Y  J N(mum2)a12fQ(cbl)fs(w2)dcbid(b2,  (D.2) 

except  that  Eq.  (D.2)  does  not  take  into  account  the  vari¬ 
able  i/f  (see  Eq.  (C.l))  and  the  fact  that  a  fragment’s  orbit 
is  often  only  partly  in  the  SOI  after  a  fragmentation  event. 
The  expected  proportion  of  time  that  a  fragment  -  of  type 
F *  and  produced  by  an  intact  of  type  a  e  U'  or  a  fragment 
of  type  Fhj  in  a  catastrophic  collision  -  spends  in  the  SOI  is 
J  °  if  t  ^  a; 

S*  l  / Wi(cu)w2(co; » c)/fl(co;  z,  k)  dm  if  t  =  a.  °'3 

The  integral  in  Eq.  (D.3)  is  the  expectation  of  the  product 
of  the  weights  wi  and  w2  defined  in  Sections  B.l  and  B.2.  A 
benign  fragment  of  type  F*  generates  only  benign  frag¬ 
ments  and  so  w2  is  excluded: 


JO  if  t  7^  a; 

\  J w\(m)fa(m-,z,b)  dm  if  z  =  a. 


(D.4) 


These  expectations  are  independent  of  what  follows  be¬ 
cause  they  concern  a  property  of  the  proportion  of  frag¬ 
ments  generated,  which,  except  for  intact  source  type,  is 
independent  of  properties  of  the  colliding  fragments. 

We  now  define  the  function  D™,  which  quantifies  -  for 
specific  fragment  property  values  -  the  number  of  frag¬ 
ments  emerging  from  a  collision  weighted  by  the  fraction 
of  time  s™  in  Eqs.  (D.3)  and  (D.4)  that  fragments  spend 
in  the  SOI.  Then,  by  Eq.  (D.2),  we  have 

K  =  VfE[Ky^\  (D-5) 


We  need  the  specific  power  law  distribution  given  in  Eq.  (4) 
in  Johnson  et  al.  (2001),  which  states  that  the  number  of 
fragments  from  a  collision  having  characteristic  length  at 
least  Lc  m  is 


N(LC;M )  =  0.1  M0J5L-ll\  (D.6) 

where  M  is  the  mass  of  the  objects  involved  if  the  collision 
is  catastrophic,  and  M  is  the  product  of  the  mass  of  the  less 
massive  object  and  the  relative  collision  velocity  vc  in  a 
noncatastrophic  collision.  We  are  interested  only  in  frag¬ 
ments  having  Lc  ^  Lc  =  1 0  cm,  and  for  convenience  we 
define  N(M)  =  N(LC:M). 

Consider  a  catastrophic  collision  between  two  objects  of 
masses  mh  i  =  1,2.  The  collision  produces  N(m\  +  m2)  frag¬ 
ments.  We  assume  -  this  assumption  goes  beyond  the 
model  in  Johnson  et  al.  (2001)  -  that  the  proportion  of 
fragments  object  i  generates  is  m  mm  ,  i  =  1,2.  Each  object 
produces  fragments  of  various  types. 


We  first  derive  F>“  for  intact-intact  collisions.  A  colli¬ 
sion  between  two  intacts  of  types  a,  y  e  U1  generates 


N(2M)  (sj 


(D.7) 


fragments  of  type  F*.  In  intact-intact  collisions,  D™  and 
are  independent. 

A  collision  between  an  intact  of  type  a.  e  V1  and  a 
benign  fragment  of  type  Fbz  and  having  mass  m  is  noncata¬ 
strophic;  only  the  benign  fragment  breaks  up.  The  frag¬ 
ment  generates 

=  N(mvc)sf  -  1  (D.8) 

fragments  also  of  type  Fb,  where  we  make  use  of 
vc  =  12.1km/s  from  Eq.  (B.8).  The  —1  accounts  for  the  loss 
of  the  original  fragment. 

A  catastrophic  collision  between  an  intact  of  type  a  e  U1 
having  mass  M  and  a  fragment  of  type  Fj  having  mass  m 
generates,  by  Eqs.  (D.2)-(D,4), 

D<, -x™)  (d.9) 


fragments  of  type  F^2,  where  <5(ti,  r2)  =  1  if  X\  =  t2  and  0 
otherwise.  The  second  term  in  parentheses  is  neglected  in 
the  fragment-fragment  version  of  the  model. 

Finally,  by  Eqs.  (D.2)-(D.4),  a  catastrophic  collision 
between  two  fragments  of  types  F*1  and  F*2  generates 

=  N(mi  +  m2)  (  — — — +  — — — sx$?  ) 

FnFn  \mi+m2  Wi  +  m2  A2  J 

—  <5(t3,  tx)  —  5(t3,  t2)  (D.10) 


fragments  of  type  F*3 . 

We  conclude  this  subsection  by  taking  the  expectation  in 
Eq.  (D.5).  For  an  intact-intact  collision,  the  expressions  do 
not  have  any  random  variables,  and 

N(2M)  Vnil/jGay 


For  a  collision  between  an  intact  and  a  benign  fragment, 
we  have 


VUEW  ^Dz^t(TxFb-  t, 

=  V  ’ 


(D.12) 


and  for  a  collision  between  an  intact  and  a  hazardous  frag¬ 
ment,  we  have 


VnEjD^  ;  Ti ,  * 

eit- — (°-13) 

Finally,  the  expectation  for  a  catastrophic  collision  be¬ 
tween  two  fragments  is 
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Table  E.l 

The  inverse  of  the  decay  rates  (in  years)  derived  in  Section  E.  The  columns  correspond  to  the  subscript;  e.g.,  ps'  =  1.10  x  10* * * 4 * * * * * 10. 

pT1  R  S  FhR  FbR  Fhs  Fbs 

1.10  x  104  1.10  x  104  1.23  x  103  3.12  x  102  2.64  x  103  7.57  x  102 


Vu^ww  5  Ti’  Kl  ’ Z2’  K2] 


(D.14) 


Although  Eqs.  (D.11)-(D.14)  yield  <5™,  we  display  values 
for  py ,  which  is  roughly  the  average  number  of  fragments 
of  a  certain  type  generated  by  a  certain  type  of  collision. 
These  numerical  values  for  both  versions  of  the  model  ap¬ 
pear  in  Tables  C.l  and  C.2. 


Appendix  E.  Decay  rates 


In  this  section,  we  estimate  the  decay  rates  for  intacts 
without  deorbit  capability  ( nR  and  jin)  and  fragments 
As  noted  in  the  main  text,  we  set  ji~l  =  3  years  for 
intacts  with  deorbit  capability  (ji0) .  We  estimate  the  decay 
rates  nR,  fi„,  and  /v  as  reciprocals  of  residence  times.  Table 

4  in  Rossi  et  al.  (1994)  gives  the  residence  time  in  years  of  a 
satellite  in  various  shells,  and  the  residence  time  for  the  SOI 

is  110  years.  The  table  is  calibrated  for  an  A/M  ratio  of 
1  m2/kg;  the  residence  time  for  a  satellite  having  a  different 

value  of  At M  is  obtained  by  dividing  the  time  by  that 

value. 

Recommended  values  for  the  A/ M  ratio  for  intacts  range 
between  2  x  10-3  and  2  x  10-2  m2/kg  (Kaula,  1983),  which 

yield  residence  times  of  between  5.5  x  103  and  5.5  x 

104  years.  Either  of  these  values  makes  the  decay  rate  for 
intacts  insignificant  in  the  next  few  centuries.  We  use 
HR  =  Hn  =  =  9.1  x  10_5/year. 

The  four  fragment  types  F *  have  different  average  A/M 
values.  By  Eqs.  (B.2)  and  (B.ll),  we  compute  the  fragment 
decay  rates  via 


[1fk  — 


Ew[10x;  r,  k] 

110 


(E.l) 


The  numerical  results  from  Eq.  (E.l)  appear  in  Table  E.l. 
Kessler  (2000)  uses  a  decay  rate  of  4g3  per  year,  which  is  be¬ 
tween  the  four  fragment  decay  rates  calculated  from  Eq. 
(E.l). 


Appendix  F.  Initial  conditions 

In  this  section,  we  allocate  the  total  number  of  initial 
fragments  Fs(0)  among  the  four  fragment  types.  In  Farinella 
and  Cordelli  (1991),  it  is  suggested  that  most  fragments  are 
from  exploding  rocket  bodies.  We  assume  90  %  of  the  initial 
fragments  are  from  rocket  bodies,  with  one  exception  that  is 
allocated  entirely  to  a  spacecraft  source:  by  separately  ana¬ 
lyzing  the  fragment  entries  in  Space  Track  (2008)  that  addi¬ 
tionally  contain  the  tag  “FENGYUN  1C”,  we  calculate  that 


400.9  of  the  955.5  initial  effective  fragments  in  the  SOI  are 
from  the  Chinese  ASAT  test  resulting  in  the  breakup  of 
the  FengYun  1C  spacecraft,  which  broke  up  at  «  850  km 
(Liou  and  Portman,  2007).  Though  the  number  90%  is  just 
a  guess,  the  results  are  insensitive  to  it;  moreover,  a  smaller 
proportion  from  rocket  bodies  yields  a  higher  fraction  of 
hazardous  fragments,  and  so  choosing  a  high  proportion 
is  conservative.  Finally,  we  use  the  parameters  s ™  in 
Eq.  (D.3),  i.e.,  the  proportion  of  fragments  of  type  Fz  gen¬ 
erated  by  the  breakup  of  an  intact  a,  to  determine  the  mix 
of  hazardous  and  benign  fragments.  This  is  a  conservative 
choice  because  benign  fragments  deorbit  faster  than  hazard¬ 
ous  ones  on  average  and  so  if  no  new  breakups  occur,  the 
proportion  of  hazardous  fragments  increases  with 
time.  These  allocations  result  in  the  initial  conditions 
Fhs( 0)  =  169.8,7^(0)  =  106.2, 7^(0)  =  286.5,  and  FbR{ 0)  = 
393.0.  Using  the  same  procedure,  we  find  that  of  the  400.9 
fragments  from  the  FengYun  1C  spacecraft,  149.1  are  haz¬ 
ardous  and  251.8  are  benign. 
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